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An isolated mass gap black hole or neutron star detected with astrometric microlensing 
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ABSTRACT 


We present the analysis of five black hole candidates identified from gravitational microlensing sur- 
veys. Hubble Space Telescope astrometric data and densely sampled lightcurves from ground-based 
microlensing surveys are fit with a single-source, single-lens microlensing model in order to measure the 
mass and luminosity of each lens and determine if it is a black hole. One of the five targets (OGLE- 
2011-BLG-0462/MOA-2011-BLG-191 or OB110462 for short) shows a significant > 1 mas coherent 
astrometric shift, little to no lens flux, and has an inferred lens mass of 1.6 - 4.4 Mo. This makes 
OB110462 the first definitive discovery of a compact object through astrometric microlensing and it 
is most likely either a neutron star or a low-mass black hole. This compact object lens is relatively 
nearby (0.70-1.92 kpc) and has a slow transverse motion of <30 km/s. OB110462 shows significant ten- 
sion between models well-fit to photometry vs. astrometry, making it currently difficult to distinguish 
between a neutron star and a black hole. Additional observations and modeling with more complex 
system geometries, such as binary sources are needed to resolve the puzzling nature of this object. For 
the remaining four candidates, the lens masses are < 2Mo@ and they are unlikely to be black holes; two 
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of the four are likely white dwarfs or neutron stars. We compare the full sample of five candidates to 
theoretical expectations on the number of black holes in the Milky Way (~ 108) and find reasonable 


agreement given the small sample size. 
1. INTRODUCTION 


Stellar-mass black holes are produced when massive 
stars collapse under their own gravity. Observations of 
black holes (BHs) are a key ingredient for understand- 
ing outstanding questions in massive stellar evolution, 
such as which stars explode, which stars produce neu- 
tron stars vs. BHs, and whether there is a gap between 
the heaviest neutron stars (NSs) and the lightest BHs. 

Black holes are abundant. There are predicted to 
be 10’ — 10° stellar-mass BHs in the Milky Way alone 
(Shapiro & Teukolsky 1983; Samland 1998; Timmes 
et al. 1996; Agol et al. 2002; Sartore & Treves 2010). 
However, only about two dozen have been definitively 
detected, all in binaries with dynamical mass measure- 
ments (Corral-Santana et al. 2016; Thompson et al. 
2019). Beyond the Milky Way, over 80 binary BH merg- 
ers have been detected via gravitational waves, with 
component masses spanning the lower mass gap ~ 3Mo 
to the lower intermediate-mass BH range ~ 100M (Ab- 
bott et al. 2020). 

These BHs are not a representative sample of the pop- 
ulation, as they are all in binary systems. While most 
massive stars exist in binary or multiple systems (Sana 
2017), the majority of the BH population is expected to 
be isolated due to the disruption of the progenitor sys- 
tems (Belczynski et al. 2004; Fender et al. 2013; Wik- 
torowicz et al. 2019). 

Isolated BHs in the Milky Way can be found and 
weighed using the technique of gravitational microlens- 
ing. When a foreground lens (e.g. BH) passes in front of 
a background source star, the source light is temporar- 
ily bent and split into two unresolved images by the lens 
mass, producing a transient photometric and astromet- 
ric signal (Paczynski 1986; Hog et al. 1995; Miyamoto 
& Yoshii 1995; Walker 1995). The characteristic cross- 
section of a microlensing event is set by the angular Ein- 
stein radius, 0g = \/«Mz(71 — 75), and depends on the 
lens mass (Mz) and the parallax of the lens (zz) and 
source (7g), where k = 4G/(1AU -c?) = 8.14 mas / Mo 
is aconstant. A ~10 Me BH in the Milky Way disk lens- 
ing a background bulge star typically has a 0g = 1-3 
milliarcseconds. 

Photometric light curves can measure the duration 
of the event, tg = Oz/[ret, where [ye is the relative 
source-lens proper motion, and the microlensing paral- 
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lax, Te = (tr — 7g)/Og. Precise, multi-epoch astrom- 
etry can measure @ directly and combined with pho- 
tometry to measure the lens mass, My = 0¢/kK7p. 

Over the past 25 years, numerous photometric mi- 
crolensing surveys have been conducted to search for 
a wide variety of lenses, including massive astrophysi- 
cal compact halo objects (MACHOs) that might make 
up dark matter, stars of all masses, and, most recently, 
exoplanets (Paczynski 1986, 1991; Griest et al. 1991; 
Mao & Paczynski 1991). Current ground-based mi- 
crolensing surveys such as the Optical Gravitational 
Lensing Experiment (OGLE, Udalski et al. 1994), Mi- 
crolensing Observations in Astrophysics (MOA, Bond 
et al. 2001), and Korea Microlensing Telescope Net- 
work (KMTNet, Kim et al. 2016) monitor hundreds of 
millions of stars toward the Galactic bulge, identify- 
ing thousands of photometric microlensing events each 
year. Photometry-only searches for stellar-mass BHs 
have been attempted (e.g. Bennett et al. (2002); Mao 
et al. (2002); Wyrzykowski et al. (2016)), but have only 
been able to identify BH candidates or place loose sta- 
tistical constraints on the BH mass function. 

In contrast to the now-routine measurements of pho- 
tometric microlensing, detections of astrometric mi- 
crolensing are still at the forefront of our technical capa- 
bilities (Lu et al. 2016; Kains et al. 2017; Rybicki et al. 
2018). Only a handful of astrometric measurements of 
the gravitational deflection of light have ever been made, 
all for nearby (< 10 pc) lenses that were astrometrically 
anticipated (Eddington 1919; Sahu et al. 2017; Zurlo 
et al. 2018) and none of which were BHs. 

If there are 10° BHs in the Milky Way, they should 
contribute only about 0.1% to the stellar mass of the 
Milky Way. In contrast, they would make up around 
1% of the Milky Way’s microlensing events due to their 
larger lensing cross section. Thus, of the thousands 
of microlensing events detected each year, a few tens 
should be due to BHs (Gould 2000; Lam et al. 2020). 
However, a 1% detection rate is akin to looking for BH 
needles in a Galactic haystack. By limiting to long dura- 
tion microlensing events with tg > 120 days, the proba- 
bility of a microlensing event being a BH rises to ~ 40% 
(Lam et al. 2020). 

We present a joint photometric and astrometric analy- 
sis of five candidate BH microlensing events. We briefly 
summarize the observations in §2 and methods in §3. A 
large non-linear astrometric microlensing signal was de- 
tected in one of the five candidates named OGLE-2011- 
BLG-0462/MOA-2011-BLG-191 (OB110462 for short) 


consistent with a compact object lens. The resulting 
mass-gap BH or NS nature of OB110462 is presented 
in §4 and discussed in a broader context in §5. Final 
conclusions are presented in §6. A complete description 
of the observations, methods, and detailed results and 
discussion for all five targets is presented in the Supple- 
mental Materials within Lam et al. (2022). 


2. OBSERVATIONS AND ANALYSIS 


Five black hole candidate microlensing events were 
first identified photometrically in the OGLE-IV (Udal- 
ski et al. 2015) and MOA (Hearnshaw et al. 2006; 
Sumi 2008) Galactic Bulge surveys. The target dis- 
cussed in detail here is OB110462, which is located at 
(17:51:40.19, -29:53:26.3) towards the Galactic bulge. 
The remaining targets properties are described in Sup- 
plemental §3.1 of Lam et al. (2022). 

The photometric light curves for each candidate span 
7 to 11 years, with approximately daily cadence except 
for seasonal gaps from November to February, and pho- 
tometric precision of ~2% for each measurement. As- 
trometric observations were obtained with the Hubble 
Space Telescope (HST) in the F606W (V-band) and 
F814W (J-band) filters. Astrometric monitoring began 
around each event’s photometric peak with a typical ca- 
dence of 1-2 times a year and an astrometric precision 
of ~ 0.3 mas per epoch. A more detailed description of 
the observations, data analysis, and multi-epoch astro- 
metric alignment is presented in Supplemental §3 and 
84 of Lam et al. (2022). 


3. MODELING METHODS 


To measure the physical properties of the lens and 
source for each event, we simultaneously fit the ground- 
based photometry and HST photometry and astrome- 
try with a point-source, point-lens (PSPL) microlensing 
model including source and lens parallax. All microlens- 
ing quantities defined in this work are in the heliocentric 
reference frame. The model free parameters describing 
the lensing geometry include tg and 7g as described in 
§1 as well as the time (to) and distance (uo) of closest 
projected approach between the source and lens in the 
heliocentric frame, and the direction of the microlensing 
parallax vector 7g. The direction of 7, is defined to 
be the same as the direction of the source-lens relative 
proper motion vector!, ie. Wz || “7. For each photo- 


1 Note that we define 1,,; to be the source-lens relative proper 
motion (Supplementary §5.1 following the convention of astro- 
metric microlensing papers (e.g. Hog et al. 1995; Miyamoto & 
Yoshii 1995; Lu et al. 2016), while the exoplanet microlensing 
community typically defines jz,..; to be lens-source relative proper 


motion. 
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metric filter, the baseline brightness (mpase) and source 
flux fraction (or blend parameter, bs rp) are fit. We also 
fit the astrometric model parameters 0 and the source’s 
parallax (7g), position at to (a@sq) and proper motion on 
the sky (wg). A complete description of the microlens 
model is presented in Supplemental 85.1 of Lam et al. 
(2022). 

The best-fit model parameters and uncertainties are 
estimated using Bayesian inference and the MultiNest 
nested sampling routine (Feroz et al. 2009; Skilling 2004, 
see Supplemental §5 of Lam et al. (2022) for details). 
There is some tension between the preferred model pa- 
rameters by the photometry vs. astrometry data for 
OB110462. When using the default weight (DW) likeli- 
hood where each data point and corresponding measure- 
ment uncertainty contributes equally to the likelihood, 
the photometry dominates as it has > 100x more data 
points than the astrometry dataset. We also fit models 
with a likelihood that give equal weight (EW) to each in- 
dependent dataset (see Supplemental §5.3 of Lam et al. 
(2022) for details). With the EW likelihood, the astrom- 
etry has much more constraining power. Results from 
both models are described in §4. 


4. RESULTS 


A large (> 1 mas) astrometric microlensing signal is 
detected in OB110462 (§4.1). The four other candidates 
are presented in Supplemental §7 of Lam et al. (2022) 
and either show no significant astrometric microlensing 
signal or have a low lens mass inconsistent with a black 
hole. 


4.1. OB110462 


From the microlensing fit, we infer that OB110462 is 
a NS or a mass gap BH, depending on the likelihood 
function adopted (see Supplemental §5.3 of Lam et al. 
(2022)). The data and model for OB110462 with the 
default weight (DW) likelihood are shown in Figures 1 
(photometry) and 2 (astrometry), and the fit posteri- 
ors are summarized in Table 1. The mass posteriors 
of the lens are shown in Figure 5. The inferred Ein- 


stein crossing time ty is 280.87*8:34 days, the microlens- 


ing parallax wp is 0.121301, the Einstein radius 0, is 
3.89t}1¢, and the lens mass Mz is 3.7973'87Mo. The 
data and model for OB110462 with the equal weight 
(EW) likelihood are shown in Figures 1 (photometry) 
and 2 (astrometry), and the fit posteriors are summa- 
rized in Table 2. The inferred Einstein crossing time 
tr is 278.56 oe days, the microlensing parallax 7 is 
0.2412, the Einstein radius 0 is 4.134)3°, and the 
lens mass My is 2.151(°'Mo. Further, we find that 
the object is located relatively nearby at 0.70-1.92 kpc 
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Figure 1. HST image (inset) and photometric lightcurve 


55800 
Time (MJD) 


(top) for OB110462. In the image, OB110462 is circled in black, 


and is shown in its unmagnified state. Observations are shown as points in blue (OGLE), red, and black (HST). The default 
weighted (DW; solid) and equal weighted (EW; dashed) maximum likelihood models are are shown along with their corresponding 
residuals (DW: middle, EW: bottom). See Supplemental §5.3 of Lam et al. (2022) for more details on the two different models. 


in the direction of the Galactic bulge and has a small 
transverse velocity of < 30 km/s. Figure 3 shows the 
on-sky lensing geometry of OB110462 inferred from the 
DW and EW likelihood models, showing the relative mo- 
tions of the source and lens with respect to each other. 

The probability that OB110462 is a dark lens is 100%, 
ruling out the possibility of a stellar lens and making 
OB110462 the first detection of a compact object with 
astrometric microlensing. Assuming there is a transition 
from white dwarfs to neutron stars at 1.2Mo and neu- 
tron stars to BHs at 2.2Mo, the relative probabilities of 
WD:NS:BH are 0:0:100 for the default weighted (DW) 
fit and 6:50:44 for the equally weighted (EW) fit. 

The microlensing fit also yields information about the 
distance and transverse velocity of the lens. The lens 
is relatively nearby at a distance of 1.47 - 1.92 kpc or 
0.70 - 1.30 kpc for the DW and EW solutions, respec- 
tively. The inferred lens velocity is < 30 km/s for both 
solutions with a slower velocity 2 - 12 km/s from the 
EW solution and a faster velocity of 21 - 27 km/s from 
the DW solution. In both cases, the velocities appear 
consistent with the compact object receiving little to no 


kick, although the line-of-sight velocity is not measur- 
able from these observations. 

Based on an analysis of the source position in the 
CMD, the OB110462 source star is around the main se- 
quence turnoff on the redder and more luminous side of 
the main sequence, suggesting it is most likely a giant or 
sub-giant star. However, a main sequence source could 
still be consistent. The relative proper motion and par- 
allax also favor a star in the near side of the Bulge (see 
Supplemental §7.6 of Lam et al. (2022) for more details). 

A point-source, point-lens (PSPL) model is not the 
end of the story for OB110462. There is no PSPL model 
which can simultaneously fit both the photometry and 
the astrometry. Specifically, the direction of [ye pre- 
ferred by the photometry and astrometry are different. 
The best fit PSPL model for the default weight (DW) 
likelihood fits the photometry very well, but leaves a 
significant ~ 0.5 mas coherent astrometric residual in 
RA (Figure 2). The best fit PSPL model for the equal 
weight (EW) likelihood leaves a significant and coherent 
~ 0.03 mag residual in the photometry, but fits the as- 
trometry in RA better than the DW likelihood model, 
although some unexplained astrometric residuals still 
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Figure 2. OB110462 astrometry, using the default weight (DW: solid) and equal weight (EW: dashed) likelihoods. Left column, 
top to bottom: RA (Aa*) vs. time with source’s unlensed motion subtracted; residuals to the maximum likelihood (MLE) model 
for Aa* vs. time fit. HST F814W astrometry data is shown in red; HST F606W astrometry data is shown in blue. The MLE 
model is shown in black. Middle column, top to bottom: Same as left column, except Dec (A6) instead of Aa*. Right panel: 
astrometry as seen on-sky, in the barycentric frame. OB110462 shows a strong > 1 mas, non-linear astrometric microlensing 


signal. 


remain. More complex microlensing geometry models, 
such as those involving a binary source or lens, should be 
explored. As mentioned in Supplemental §4.2.5 of Lam 
et al. (2022), we apply a constant positional offset to 
the F606W data in order to make it match up with the 
F814W data. However, this filter dependent positional 
difference may actually be astrophysical and consistent 
with a small contribution from a faint companion to the 
source. Either way, both solutions indicate a NS or BH 
detection. 

The alternative explanation to the tension between 
the photometry and astrometry of OB110462 is some 
type of systematic error in one or both sets observations. 
This possibility is discussed in Supplemental §8.3 of Lam 
et al. (2022). 


4.2. Is OB110462 a BH or a NS? 


One means of further understanding the BH or NS 
nature of OB110462 would be to detect electromagnetic 
radiation from the lens. We searched existing X-ray and 
pulsar catalogs at the position of OB110462 and did 
not find any counterpart (see Supplemental §8.2 of Lam 
et al. (2022) for details). Unfortunately, OB110462 is 
not in the Gaia EDR3 catalog as it is too faint. 

Future observations of OB110462 will be useful to de- 
termine its true nature. Continued astrometric mon- 
itoring, including the remaining data from HST Cycle 
29 program GO-16760 (Lam & Lu 2021) will continue to 
improve the lens mass estimate as described in Supple- 
mental §8.2 of Lam et al. (2022). High contrast imaging 
observations in the next 5-10 years would also be worth- 
while. If the lens is indeed a solitary NS or BH, lack of 
an optical/infrared lens detection would bolster support 
of the dark isolated lens interpretation, but be unable 
to distinguish between NS or BH as the relative proper 
motion differentiating the DW and EW solutions cannot 
be measured in the case of a dark lens. In that case, only 
a deep, targeted X-ray observation could help in differ- 
entiating between the NS vs. BH scenarios. On the 
other hand, optical/infrared detection of a lens separate 
from the source would point to a binary lens scenario. 
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Figure 3. On-sky lensing geometry of OB110462 showing the resolved motion of the lens and source, as inferred from the 
default weight (DW) model (left) and the equal weight (EW) model (right); see Supplemental §5.3 of Lam et al. (2022) for 
details about the different models. The Einstein ring is shown as a gray circle of radius 0g. The solid red line shows the 
trajectory of the source, while the dashed black line shows the trajectory of the lens. Note that the red line shows the unlensed 
position or the source, and not the centroid of the source’s lensed images. The dots on top of the trajectories are spaced at 
intervals of 100 days. The red and black arrows indicate the proper motion of the source and lens, respectively. The tail of the 
arrow is at the location of the source and lens at time to; the length of the arrows are proportional to the magnitude of the 


source and lens proper motions. 


Binary source or binary lens models can be further ex- 
plored through precision radial velocity searches. 


4.3. Number of detected BHs 


Next, we compare our observed BH yield to the the- 
oretical expectation calculated using the PopSyCLE sim- 
ulations. PopSyCLE’s Galactic model contains 2 x 108 
BHs ranging from 5 — 16Mo (Lam et al. 2020). For a 
sample of simulated events that would be observable by 
OGLE (see Table 4 of Lam et al. (2020)), we calculate 
the fraction of those events due to BHs as a function 
of the Einstein crossing time tg as described in Supple- 
mental §7.9 in Lam et al. (2022). Figure 6 shows the ex- 
pectation of detecting Ney = 0,1,...,5 BHs within our 
sample of 5 targets. The probability of detecting 0 or 1 
BHs in our simulation is ~25% and ~45%, respectively. 
This estimate is consistent with our single detection of 
a NS-BH object. Note that in PopSyCLE, there are no 
2— 5M. mass gap NSs or BHs in the simulation, and 
hence no exact OB110462 analogue. 


5. DISCUSSION 


OB110462 is the first definitive detection of a com- 
pact object discovered with astrometric microlensing. 
Depending on the likelihoood function used to evalu- 
ate the fit (see Supplemental §5.3 of Lam et al. (2022)), 


it is either a neutron star (50% probability for the EW 
likelihood), a black hole (44% probability for the EW 
likelihood, 100% probability for the DW likelihood), or 
a white dwarf (6% probability for the EW likelihood). 
The other four candidates discussed in detail in Supple- 
mental §7 of Lam et al. (2022) are mostly likely stars, 
white dwarfs or neutron stars, although a NS-BH mass 
gap object cannot be ruled out in two cases. 

Here we discuss the observed BH yield as compared 
to theoretical expectations (§5.1) and and implications 
for the BH mass function (85.2). 


5.1. Comparison to simulations 
Oiled. Te-tE-Oc,max 


As described in Lam et al. (2020), BH candidates 
can be identified photometrically by their long tg and 
small tz, and confirmed astrometrically by measuring 
the maximum astrometric shift dcmax = On/\/8 (see 
Supplemental §5.1 ). Figure 4 shows the 1 — 2 — 30 
posterior contours of wg vs. ty (left) and mH Vs. bc,max 
(right) of the microlensing models for the 5 targets, com- 
pared against simulated microlensing events generated 
by the PopSyCLE software (Lam et al. 2020). By com- 
paring the mg — tg and 7g — b¢,maz posteriors against 
simulation, we can gain a more intuitive understanding 
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Figure 4. Microlensing parallax 7g vs. Einstein crossing time tm (left) and maximum astrometric shift dcmaz (right). Points 
are from the PopSyCLE simulation. Contours are 1 — 2 — 30 (39.3-86.5-98.9%) credible regions from the microlensing model fits 
to the five BH candidates. There are two fits for OB110462 (default weight (DW) and equal weight (EW); see Supplemental 
§5.3 of Lam et al. (2022) for details). The OB110462 DW solution has a smaller mz than the OB110462 EW solution, and 
has a correspondingly more massive lens mass. Both solutions fall solidly in the NS-BH mass gap, making OB110462 the best 
BH-candidate. MB09260 and OB110310 are most likely white dwarfs or neutron stars, although due to uncertainty in mm and 
6c,max higher and lower mass lenses cannot be definitively ruled out. OB110037 and MB10364 are not BHs as they have very 


large 7p, as well as relatively short tg and small d¢,maz. MBO09260, 


in the Supplemental paper. 


of the inferred lens types for the targets (Supplemental 
Table 9 of Lam et al. (2022)). 

Both EW and DW models for OB110462 fall solidly 
within the 2— 5M. mass gap shown in the tg — d¢,mazx 
parameter space. Because the EW solution leads to a 
larger and more uncertain value of tz than the DW 
solution, a neutron star or even white dwarf lens is a 
possibility. On the other hand, the DW solution prefers 
a smaller and more well constrained value of 7g than 
the EW solution, leading to a much more definitive so- 
lution of a mass-gap BH. However, both the EW and 
DW solution for OB110462 fall in a somewhat unusual 
part of the 7_—tg parameter space for BHs: typical BH 
7m are around 0.02, while for OB110462, mz is around 
0.1. This is because PopSyCLE simulations only contain 
BHs with masses from ~ 5—16Mo; if OB110462 is truly 
a mass-gap BH, it would not correspond to any BHs in 
the simulation. 


5.2. OB110462 in comparison to the BH population 


Several attempts have been made to determine the 
Milky Way BH mass function using dynamical mass 
measurements of BHs in binaries. The mass-gap be- 
tween 3-5 Mo was first observed in low-mass X-ray bi- 


MB10364, OB110037, and OB110310 are discussed in detail 


naries (Bailyn et al. 1998; Ozel et al. 2010; Farr et al. 
2011). However, more recent detections of BHs in this 
mass range from both gravitational wave mergers and 
in non-interacting binaries suggests that the mass-gap 
may actually be filled with BHs (e.g. Abbott et al. 2017; 
Thompson et al. 2019). For a complete description of 
the Milky Way BHs found to date, see Supplemental 
§8.1 of Lam et al. (2022). If OB110462 is a BH, it too 
shows that the mass function of BHs extends into the 
mass-gap regime. 

It is somewhat surprising that massive BHs have not 
been found in our microlensing search as well as in 
searches for BHs in wide binaries (Rowan et al. 2021; 
Thompson et al. 2019; El-Badry et al. 2022). Both mi- 
crolensing and radial velocity searches should be biased 
toward finding 10 Me objects more easily than 3 Mo 
objects as described in Supplemental §8.1 of Lam et al. 
(2022). It may be that the selection bias is offset by 
the steep mass function for massive stars and thus BHs. 
As the sample of BHs in the Milky Way grows, a more 
quantitative analysis of the sample selection will be es- 
sential to constrain the true BH mass function. 


6. CONCLUSION 


Table 1. OB110462 DW Fit Values 
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Mbase,o (mag) 16.4170 :0nt 16.41 16.41 
bsFr,HE8 aso"? 0.89 0.91 
Mbase,Hs (mag)  19.8679-0°6 19.86 19.87 
bsFF,H6 aga 0.95 0.94 
Mbase,He (mag) 22.0479 O08 22.05 22.04 
Mz (Mo) a1 3.01 3.58 
mr, (mas) apo ee 0.48 0.67 
Tre (mas) 0.4813-03 0.36 0.55 
it,B (mas/yr) —2.6410-38 -2.54 -2.28 
t.n (mas/yr) Lage 0.26 1.69 
Lret,p (mas/yr) g.aotes 0.28 0.03 
Lret,w (mas/yr)  —5.0243-% -3.81 -5.26 
0% (mas) 3.89512 2.98 4.00 
TE Ge 0.14 0.13 
Jc,max (mas) aaa 1.05 1.41 








NotE—The columns list the median +1lo (68%) credible 
intervals, maximum a posteriori (MAP) solution, and 
and maximum likelihood estimator (MLE) solution. 





We analyze five microlensing events with candidate 
BH lenses. Combining HST astrometry and densely 
sampled ground-based photometry, we derive masses for 
these five lenses as well as their probability of being a 
BH. Of the five targets, we make one definitive > 1 mas 
detection of astrometric microlensing (OB110462). The 
mass of the lens in OB110462 is in the range 1.6-4.4 
Mo, making it the first detection of a compact object 
through astrometric microlensing. 

We use our detection of a mass-gap BH or neutron star 
and the non-detections of BHs in the rest of the sam- 
ple to observationally constrain the number of BHs in 
the Milky Way. Our observational BH yield currently 
agrees with simulations assuming 2 x 10° BHs in the 
Milky Way, albeit with very large uncertainties due to 


Table 2. OB110462 EW Fit Values 
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the small sample size. The ability to place more strin- 
gent constraints on the number and mass distribution of 
Galactic BHs will require larger samples, such as those 
that may be delivered by the Roman Space Telescope’s 
microlensing survey. 

Astrometric microlensing holds the key to uncovering 
the hidden BH population. Further pursuit and refine- 
ment of the event selection, observing, and modeling 
process will fulfill the full promise of this technique and 
its ability to reveal the properties of Galactic BHs. 


Shortly prior to this work being submitted for review, 
we learned of an independent analysis of OB110462 car- 
ried out by Sahu et al. (2022). Notably, they reach a 
different conclusion about the mass of the lens (7.1 + 
1.3Mo). It is not clear whether the discrepancy is due 
to the use of different datasets (e.g. we include an addi- 
tional epoch of 2021 HST data), performing the analysis 
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Figure 5. Lens mass posterior probabilities for the five mi- 
crolensing BH candidates. Two mass posteriors are shown 
for OB110462, one for each model (default weight (DW) and 
equal weight (EW)). See Supplemental §5.3 of Lam et al. 
(2022) for details on the two models. MB09260, MB10364, 
OB110037, and OB110310 are discussed in detail in the Sup- 
plemental paper. 
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differently (e.g. we explore solutions allowed by both 
photometry and astrometry using different likelihood 
weights), or a combination of both. In addition, al- 
though both analyses make clear detections of an as- 
trometric deflection, the direction of the deflections are 
in opposing directions in RA. Preliminary work shows 
that different choice of reference stars across the two 
teams is not the source of the discrepancy. However, 
significant further work is required to fully understand 
the differences between the two analyses. 
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ABSTRACT 


This Supplement provides supporting material for Lam et al. (2022). We briefly summarize past 
gravitational microlensing searches for black holes and present details of the observations, analysis, and 
modeling of five black hole (BH) candidates observed with both ground-based photometric microlensing 
surveys and Hubble Space Telescope astrometry and photometry. We present detailed results for four 
of the five candidates that show no or low probability for the lens to be a black hole. In these cases, 
the lens masses are < 2Mo and two of the four are likely white dwarfs or neutron stars. We also 
present detailed methods for comparing the full sample of five candidates to theoretical expectations 
of the number of BHs in the Milky Way (~ 108). 

1. SUPPLEMENTAL INTRODUCTION This paper is a supplement to Lam et al. (2022), in 
which we present results from a search for stellar mass 
black holes in the Milky Way using gravitational mi- 


Corresponding author: Casey Y. Lam crolensing. Lam et al. (2022) and this supplement de- 
casey_lam@berkeley.edu scribe our analysis of five Milky Way black hole (BH) 
candidates, which constitute a sufficient sample to place 


* OGLE collaboration 
+ MOA collaboration preliminary constraints on the number of isolated black 


holes in the Milky Way. 


This supplement is organized as follows. In §2, we 
present a detailed review of past microlensing searches 
for BHs. In §3, the datasets and reduction processes for 
our BH search are described, and in §4 the photometric 
and astrometric analysis of the HST data is explained. 
In §5, the microlensing modeling and fitting procedure 
are described in detail, while §6 describes how to com- 
bine the high-resolution imaging and results of the mi- 
crolensing modeling to constrain the lens’ luminosity. §7 
presents the results of the modeling for the five candi- 
dates. Note that results for the single mass-gap black 
hole or neutron star candidate are presented in the main 
paper Lam et al. (2022). Details on constraints on the 
Milky Way BH population and future BH microlensing 
searches are discussed in 88. 


2. PAST MICROLENSING SEARCHES FOR 
BLACK HOLES 


The advent of ground-based microlensing surveys 
provided a new avenue to search for isolated black 
holes (BHs). Notable BH candidates identified with 
photometric microlensing include MACHO-96-BLG-5, 
MACHO-98-BLG-6 (Bennett et al. 2002), and MACHO- 
99-BLG-22/OGLE-1999-BUL-32 (Mao et al. 2002). 
Poindexter et al. (2005) found that MACHO-99-BLG- 
22 is likely a BH, MACHO-96-BLG-5 is possibly a 
BH, and MACHO-98-BLG-6 is most likely not a BH, 
with BH lens probabilities of 78%, 37%, and 2%, re- 
spectively. On the other hand, Abdurrahman et al. 
(2021) found that MACHO-96-BLG-5 and MACHO-98- 
BLG-6 are still good BH candidates, ruling out non-BH 
lenses for source-lens relative proper motions larger than 
2.5 mas/yr. However, mass estimates for these lens- 
ing events cannot be made without invoking a Galac- 
tic model. These candidates only had photometric mi- 
crolensing observations, which alone cannot constrain 
the mass of the lens, unless rare higher-order effects such 
as finite-source effects are detected. 

As mentioned briefly at the end of Paczynski (1986), 
microlensing also has an astrometric signature, in which 
the centroid of the image is displaced from the source’s 
true position (Hog et al. 1995; Miyamoto & Yoshii 1995; 
Walker 1995). In contrast to the now-routine measure- 
ments of photometric microlensing, detections of astro- 
metric microlensing are still at the forefront of our tech- 
nical capabilities. Typical astrometric shifts toward the 
Bulge are O(0.01—1) mas, and few existing facilities are 
currently capable of the astrometric precision to perform 
this measurement. Only a handful of astrometric mea- 
surements of the gravitational deflection of light have 
ever been made, all for nearby (< 10 pc) lenses that 
were astrometrically anticipated (Eddington 1919; Sahu 


et al. 2017; Zurlo et al. 2018). However, a combination of 
photometric and astrometric microlensing together can 
determine the mass of the lensing object, making de- 
tection of astrometric microlensing important for BH 
searches. 

To date, there have been two endeavors to measure 
lens masses by combining photometric and astrometric 
microlensing. Lu et al. (2016) attempted a measure- 
ment with Keck laser guide star adaptive optics (LGS 
AO), but no detections of astrometric microlensing were 
made. Kains et al. (2017) reported a detection of astro- 
metric microlensing made with the Hubble Space Tele- 
scope (HST), but the signal was very weak and no lens 
masses were well constrained. The currently operating 
Gaia mission is also anticipated to make measurements 
of astrometric microlensing, which can be searched for 
once per-epoch astrometry is released (McGill et al. 
(2020) and references therein). A handful of these de- 
flections should be due to BHs (Rybicki et al. 2018). 


3. OBSERVATIONS IN DETAIL 
3.1. Targets and Selection Criteria 


Five BH candidate microlensing events! were selected 
from the OGLE Early Warning System? (Udalski et al. 
2015) and MOA Alerts? to be imaged with HST. These 
targets were selected to have long (tz > 200 days) dura- 
tion, no light contribution from the lens, and high mag- 
nifications to allow detection of parallax signals, making 
them good isolated BH candidates (Sahu 2009). Prelim- 
inary results by Sahu et al. (2017) reported that all five 
candidates were low-mass (< 0.5M,) stars, with no BH 
detections in the sample. We chose to re-analyze these 
targets in order to use both BH detections and non- 
detections in the sample to constrain the total number 
of BHs in the Milky Way. 

Three of the targets, OGLE-2011-BLG-0037/MOA- 
2011-BLG-039, OGLE-2011-0310/MOA-2011-BLG-332, 
and OGLE-2011-BLG-0462/MOA-2011-BLG-191 (here- 
after OB110037, OB110310, and OB110462), were 
alerted by both OGLE and MOA. The other two targets, 
MOA-2009-BLG-260 and MOA-2010-BLG-364 (here- 
after MB09260 and MB10364), were only alerted by 
MOA. Table 1 lists their coordinates. Figures 3 - 8 show 
the lightcurves of the targets. 


3.2. MOA 


1 Three additional other targets were initially observed with HST, 


but dropped from the target list after a year (Sahu 2012). 
? http: //ogle.astrouw.edu.pl/ogle4/ews/ews. html 


3 https: //www.massey.ac.nz/~iabond/moa/alerts/ 


Table 1. Target Summary 











Short Name OGLE Alert Name MOA Alert Name RA (J2000.0) Dec (J2000.0) 
MB09260 i MOA-2009-BLG-260 = 17:58:28.561 -26:50:20.88 
MB10364 = MOA-2010-BLG-364 17:57:05.401 -34:27:05.01 
OB110037 OGLE-2011-BLG-0037 MOA-2011-BLG-039 — 17:55:55.83 -30:33:39.7 
OB110310 OGLE-2011-BLG-0310 MOA-2011-BLG-332  17:51:25.39 -30:24:35.0 
OB110462 OGLE-2011-BLG-0462 MOA-2011-BLG-191  17:51:40.19 -29:53:26.3 





The MOA-II survey is carried out with a 1.8m 
telescope at Mt. John University Observatory in New 
Zealand (Hearnshaw et al. 2006; Sumi 2008). The see- 
ing at the site ranges from ~ 1.8” — 3.5”, with the me- 
dian seeing being ~ 2.5’. The telescope has a 2.2 deg? 
FOV, with a 10-chip CCD camera with a plate scale 
0.57” /pixel. The main observations are taken using the 
MOA-Red (630-1000 nm) filter (Bond et al. 2001). 

MOA data for MB09260 and MB10364 were re- 
duced as described in Bond et al. (2017).4. The MOA 
lightcurves are photometrically calibrated to the OGLE- 
ITI I-band. 


3.3. OGLE 


The OGLE-IV survey is carried out at the 1.3-m War- 
saw telescope at Las Campanas Observatory in Chile 
(Udalski et al. 2015). The seeing at the site ranges from 
~ 1.0” — 2.0’, with the median seeing being ~ 1.3”. 
The telescope has a 1.4 deg? FOV, with a 32-chip CCD 
camera with a plate scale of 0.26” /pixel. The main ob- 
servations are taken using the OGLE-I filter, which is 
similar to Cousins L-band. The data was reduced us- 
ing the Difference Image Analysis technique as imple- 
mented by Wozniak (2000). In addition, we rescaled the 
photometric uncertainties of OB110462 according to the 
method described in Skowron et al. (2016). 

OGLE data is only available for OB110037, 
OB110310, and OB110462. The magnification of 
MB09260 was not observed by OGLE as it occurred 
during the OGLE-III to OGLE-IV upgrade. MB10364 
is located in a gap in the detectors of the OGLE camera. 


3.4. HST 


HST observations come from a multi-year campaign 
following up these five targets (GO-11707, GO-12322, 
GO-12670, GO-12986, GO-13458, GO-14783; PI: K. C. 
Sahu). Observations were taken with the UVIS chan- 
nel on the Wide Field Camera 3 (WFC3) in two dif- 


4 0B110037, OB110310, and OB110462 also have MOA 
lightcurves. For simplicity we only present the OGLE lightcurve 
fits for those events, since the seeing at OGLE is better. 


ferent wide-band filters, F606W (V-band) and F814W 
(I-band). Table 2 summarizes the HST observations. 

The WFC3 UVIS channel is composed of two 2kx 4k 
CCDs and has a 162” x 162” field of view with a plate 
scale 0.04” /pixel. WFC3 UVIS supports sub-arraying, 
in which only a portion of the entire detector is read out, 
which can reduce data volume or exposure time and in- 
crease observational efficiency. All observations prior to 
2011-07-22 were taken with the UVIS1-2K4-SUB subar- 
ray mode. Beginning HST Cycle 18, more subarray sizes 
were made available, and observations after 2011-07-22 
were taken with the UVIS2-2K2C-SUB subarray mode, 
a 2kx2k subarray. 

Additional observations of OB110462 commenced in 
Cycle 29 (GO-16760; PI: C. Lam). These were taken 
in as similar a configuration as possible to the later 
epochs of the archival program, using WFC3 UVIS in 
UVIS2-2K2C-SUB subarraying mode, with observations 
in F606W and F814W filters. The first set of observa- 
tions from this program was taken October 2021 and 
is presented here; an additional set of observations is 
anticipated to be taken Fall 2022 (Lam & Lu 2021). 


3.5. Gaia 


Gaia is an all-sky scanning astrometric space mission 
(Gaia Collaboration et al. 2016). All of the targets, with 
the exception of OB110462, are found in Gaia Early 
Data Release 3 (EDR3, Table 4). Gaia EDR3 covers 
the period from 25 July 2014 to 28 May 2017 (Gaia 
Collaboration et al. 2021). MB10364, OB110037, and 
OB110310 have proper motions and parallaxes, while 
MB09260 does not. OB110462 was not in Gaia as it is 
too faint. We note that there is a Gaia source located 
~ 0.35” away from OB110462 (Gaia EDR3 Source ID 
4056442477683080960), which coincides with the bright 
star directly west of the target seen in the HST images 
(Figure 1). 


4. HST DATA ANALYSIS 
4.1. Reduction 
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Central 3” x 3” of HST WFC3-UVIS F814W combined images of the observed fields, centered on the target 


(circled). These images are of the target at or near baseline, i.e. unmagnified. The color stretch is logarithmic. Note that the 


color scale is not the same across panels. 
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Figure 2. CMDs for each field. The target at baseline magnitude and color is marked as a red star. 


The HST archival data was accessed from the Mikul- 
ski Archive for Space Telescopes? (MAST) in June 2021. 
The Cycle 29 data was accessed October 2021. For 
the following analysis, we employed the calibrated, flat- 
fielded, individual exposures corrected for charge trans- 


5 https://archive.stsci.edu/hst / 


fer efficiency (CTE, HST files with suffix _f1c°). The 
archival data was processed with Version 3.6.0 (Dec-31- 
2020) of the calwf3 pipeline, using Version 2.0 of the 
CTE correction algorithm. 


6 See Gennaro (2018) for a full description of the different file name 


suffixes. 
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Figure 3. Top panel: Detrended MB09260 lightcurve, as seen by MOA and HST F814W and F606W. The maximum likelihood 
model (MLE, described in §7) is plotted over the data. Second from top panel: The residuals to the MLE model. The Gaussian 
Process (GP) model is plotted on top of the residual. We emphasize the residual is not independently fit by the GP, but is 
simultaneously fit with the model parameters; this is purely to visualize the data (also see Golovich et al. (2022)). See §5 
for more details about the fitting procedure. Second from bottom panel: Same as top panel, but zoomed into the three most 
magnified years. Bottom panel: Same as second from top panel, but zoomed into the three most magnified years. 


CTE can alter astrometry at the milliarcsecond level, 
hence it is important to use _f1c files. However, even the 
_flc files do not necessarily fix all problems associated 
with CTE (Kuhn & Anderson 2021). Exploration of 
other methods of CTE correction will be explored in 
future work. At the present we mitigate CTE effects via 
other methods (§4.2.2) and validate our astrometry to 
ensure it is not distorted by CTE. 

Images were converted into calibrated star lists via the 
following steps. 


1. Star list extraction from individual frames. Star 
lists were extracted from the individual -f1lc ex- 
posures by modeling the PSFs of sources with 
hstipass, an updated version of the software de- 
scribed in Anderson & King (2006). Empirical 
filter-dependent PSF models as described in An- 
derson (2016) and geometric distortion solutions 
as described in Bellini et al. (2011) were used when 
performing source extraction with hstipass. 


2. Combined star list for one epoch. Within a single 
epoch and filter, multiple star lists were aligned to 


a common coordinate system in an iterative man- 
ner using xym2mat and xym2bar (Anderson & King 
2006), which include the distortion solution for the 
WFC3 camera and filters (Bellini et al. 2011) to 
produce a single matched star list. 


3. Photometric calibration. Lastly, a zero-point is ap- 
plied to the star lists to convert from instrumen- 
tal to Vega magnitudes. Star lists were calibrated 
against photometrically calibrated star lists on the 
Hubble Legacy Archive, Data Release 10 (HLA 
DR10). A magnitude offset is applied later during 
the astrometric alignment (§4.3) to obtain more 
precise relative photometry. 


Note that data taken in F606W and F814W filters are 
treated as independent measurements. That is, obser- 
vations taken on the same date are treated as distinct 
epochs, and are not combined into a single star list, as 
the importance of filter dependence in astrometry is not 
well established. See §4.2.5 and Appendix H for further 
details. 
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Figure 4. Same as Figure 3, but for MB10364. Note we do not include GP in the MB10364 fit. 
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Figure 5. Same as Figure 3, but for OB110037. Instead of MOA data, the blue data are OGLE data. 


Certain epochs were excluded from the analysis; these e Epochs with only a single frame. Observations 
are marked with an asterisk in Table 2. The reason for with only a single frame per filter cannot pro- 
their exclusion is detailed as follows. 
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Figure 6. Same as Figure 5 


duce any useful photometric or astrometric con- 
straints using hstipass. This is the case for 
the MB09260 2009-10-01 F814W and 2009-10-19, 
2010-03-22, and 2010-06-14 F606W epochs. In the 
MB09260 2010-10-20 F606W epoch, a cosmic ray 
in one of the exposures interfered with the extrac- 
tion of the target, effectively leaving only a single 
usable frame. 


Multiple exposure times. Although in principle 
mixing multiple exposure times in a single epoch 
is possible, in practice most of the data was ob- 
tained with several long exposures and only a sin- 
gle frame with a shorter exposure. Rather than 
analyze the few short frames with different detec- 
tion thresholds, PSF reference stars, and astro- 
metric reference stars, which can lead to system- 
atic errors, we choose to only use frames with the 
same exposure times within an epoch. For this 
reason, some frames from the MB09260 2009-10-01 
F606W, MB10364 2010-09-13 F814W, OB110037 
2011-08-15 F814W, OB110310 2013-10-21 F814W, 
and OB110462 2011-08-08 F606W and F814W 
epochs were not used. 


Saturation of target. No useful astrometric limits 
can be placed when the target is saturated. The 





Time (MJD) 


, but for OB110310. 


target is saturated in both filters in the MB10364 
2011-04-13 and 2011-07-22 epochs. 


Telescope pointing issues. The observations of 
OB110462 on 2017-08-11 suffered a telescope drift- 
ing issue, resulting in streaked images. 


Astrometric alignment systematics. Although 
there are no standalone issues with the observa- 
tions of OB110462 on 2013-05-13, astrometric sys- 
tematics are apparent in the reference stars when 
this epoch is astrometrically aligned along with 
the other epochs using the methodology described 
in §4.2. This is due to the difference in position 
angle of the observations taken, as the 2013-05-13 
epoch was taken at PA = 99.9 deg, while all the 
other epochs were taken with a PA different by 
~ 180 deg, with PA = 255.2 - 276.1 deg. Thus, 
the 2013-05-13 epoch is left out of the analysis. 


The other targets (MB09260, OB110037, 
OB110310) with ~ 180 deg differences in PA 
across observations do not suffer this same prob- 
lem, as there are multiple observations at each PA. 
This allows the systematics due to the ~ 180 deg 
PA flip to be calibrated out during the astrometric 
alignment. 
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Figure 7. Same as Figure 5, but for OB110462. This maximum likelihood model was calculated using the default weighted 


likelihood described in §5.3. 
4.2. HST Astrometric Analysis 


The positional measurements extracted from the dif- 
ferent epochs of HST data (§4.1) must be transformed 
into a common reference frame in order to derive the mo- 
tion of the target. This is an iterative process, with mul- 
tiple “passes” at refining the reference frame, allowing 
for the best relative astrometry possible to be extracted. 
We follow a similar procedure as described in §4.2 of Lu 
et al. (2016). The final positions and magnitudes of the 
targets resulting from the analysis in this section is pre- 
sented in Table 7. The photometry is shown in Figures 
3 - 8 and the astrometry is shown in Figures 14 - 19. 


4.2.1. Alignment procedure 


Following standard image processing techniques, a 2- 
D polynomial transformation of the form 


(1) 
(2) 


x’ =ag tax + agy + agx? + agry + asy? t+... 





y’ = bo t byx t boy t b3x? { baxy t bey? +... 


is applied to the images in order to match them to a 
reference image. A first order 2-D polynomial transfor- 
mation (2’ = agp + a,2 + agy; y’ = bo + biz + doy) is 


an affine transformation’, which can be used to model 
translation, rotation, scaling, and shearing introduced 
by the camera. A higher-order polynomial can correct 
for additional distortions, but going beyond second or- 
der generally does not improve results as the number of 
free parameters quickly increases and results in overfit- 
ting. 

In the first pass, the HST images are aligned to the 
absolute reference frame of Gaia with a first order 2- 
D polynomial transformation to roughly establish the 
transformation. The Gaia EDR3 catalog is matched to 
the HST catalog using the pattern matching algorithm 
of Groth (1986). 

In subsequent passes, the HST images are aligned to 
themselves, using a 2-D polynomial transformation go- 
ing up to second order. It is empirically determined that 
3-4 passes gives optimal results. In each successive pass, 
the HST images are aligned to the reference frame de- 
rived in the previous pass, which continually refines the 
reference frame and derived proper motions. 


7 An affine transformation maps points to points, lines to lines, 
planes to planes, and so on. Affine tranformations preserve 
collinearity and ratios of distances. Parallel lines also remain 
parallel after an affine transformation. 
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Figure 8. Same as Figure 5, but for OB110462. This maximum likelihood model was calculated using the equal weight 


likelihood described in §5.3. 


To calculate the optimal transformation, a set of 
reference stars r& (to) are selected from the stars in 
reference frame R(to) observed at time tp. The refer- 
ence stars r® (to) are matched to corresponding stars 
uje(t) in the untransformed U(t) frame observed at 
time t. The transformation T : U(t) > R(t) is found 
by least-squares minimization of the x and y position 
residuals from the alignment 


Tres = 2 We,i Ge (t) ~~ Tee ke) (3) 


Yres = a Wy,iYma(t) — T(ywa(t)))? (4) 


where wW(,y),; is the weight for the i-th reference star, 
(xP (t),a,y®()ra) and (Y(t), y¥,(t)) are the positions 
of the stars in the reference R(t) and untransformed U(t) 
frames at time t. The positions of the reference stars 
re y(t) in the reference frame R at time ¢ are propagated 
from time to using the proper motions 


where tee are the proper motions of the i-th ref- 
erence stars in reference frame R. After the trans- 
formation T is derived, it is applied to all the stars 
u” in the U frame to obtain a transformed star list 
u®(t) = T(u"(t)), where the stars u(t) are now in 
the frame R. This yields starlists for all N epochs 
uf(t,), ...,u%(ty), where the positions of all the stars 
are now in the same reference frame R. For each star j, a 
proper motion is derived by finding the best-fit straight 
line via non-linear least squares through the n < N ob- 
servations.® 

Lu et al. (2016) examined several different weighting 
schemes and showed that the resulting astrometry is not 
affected. We choose to use weights w(z,y) i = 1/0(a,y),i5 
where o(,y),; are the positional uncertainties of the stars 
in the untransformed frame. For the positional uncer- 
tainties, instead of using the root mean square values 
Orms directly returned by hstipass, we follow Hosek 
et al. (2015) and use the error on the mean ogus/VN 
where WN is the number of frames the source is detected 
in, with an additional empirical additive error. The un- 


8 Some stars are not detected in all epochs, which is why it is 
possible to have n < N observations. 
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Table 2. HST WFC3-UVIS Observations 





























Target Epoch (UT) PA Filter Teap Nim 
(yyyy-mm-dd) (deg) (sec) 

MB09260 2009-10-01 275.0 F606W _ 10.0 4* 
100.0 2 

F814W 20.0 7 

2009-10-19 275.0 F606W 310.0 1* 
F814W 72.0 6 

2010-03-22 95.0 F606W 300.0 1* 
F814W 280.0 5 

2010-06-14 95.0 F606W 200.0 1* 
F814W 275.0 5 

2010-10-20 270.0 F606W 275.0 2* 
F814W 275.0 4 
2011-04-19 90.0 F606W 275.0 2 
F814W 275.0 4 
2011-10-24 270.0 F606W 275.0 2 
F814W 275.0 4 
2012-09-25 270.0 F606W 275.0 2 
F814W 275.0 4 
2013-06-17 105.5 F606W 275.0 2 
F814W 275.0 4 
MB10364 2010-09-13 270.0 F606W _ 1.0 1 
2.0 1 
F814W — 3.0 5 

1.0 * 
2010-10-26 277.4 F606W~ 12.0 2 
F814W_ 12.0 6 

2011-04-13 90.0 F606W 260.0 2* 

F814w 120.0 4* 

2011-07-22 260.5 F606W 160.0 4* 

F81i4w 90.0 4* 
2011-10-31 278.8 F606W — 30.0 5 
F814W — 30.0 6 
2012-09-25 270.9 F606W — 30.0 5 
F814W — 30.0 6 
2013-10-24 277.0 F606W 40.0 5 
F814W 40.0 6 
OB110037 2011-08-15 245.0 F606W — 30.0 4 
F814W 40.0 4 

20.0 3" 
2011-09-26 270.8 F606W — 30.0 5 
F814W 20.0 6 
2011-11-01 276.1 F606W 50.0 5 
F814W — 30.0 5 
2012-05-07 98.1 F606W _ 80.0 4 
F814W _ 60.0 5 
2012-09-25 270.8 F606W_ 80.0 4 
F814W _ 60.0 5 
2013-10-21 274.3 F606W 80.0 4 
F814W _ 60.0 6 
2014-10-26 275.1 F606W_ 60.0 4 
F814W 55.0 6 
2017-03-13 90.0 F606W 60.0 3 
F814W 55.0 6 
2017-09-04 256.9 F606W_ 60.0 3 
F814W 55.0 6 





NotEe— Asterisk (*) denotes observations excluded from analysis. 


Table 3. HST WFC3-UVIS Observations 














Target Epoch (UT) PA Filter Teap Nim 
(yyyy-mm-dd) (deg) (sec) 
OB110310 2011-09-21 270.0 F606W ~ 75.0 4 
F814W 75.0 5 
2011-10-31 276.5 F606W 280.0 3 
F814W 200.0 4 
2012-04-24 96.0 F606W 280.0 3 
F814W 230.0 4 
2012-09-24 271.3 F606W 280.0 3 
F814W 230.0 4 
2013-10-21 274.8 F606W 280.0 3 
F814W_ 68.0 1* 
230.0 4 
2017-03-14 90.4 F606W 270.0 3 
F814W 230.0 4 
2017-09-01 268.4 F606W 270.0 3 
F814W 230.0 4 
OB110462 2011-08-08 270.0 F606W_ 60.0 1* 
75.0 3 
F8l4w 120.0 1* 
60.0 ies 
75.0 3 
2011-10-31 276.1 F606W 280.0 3 
F814W 200.0 4 
2012-09-09 269.5 F606W 290.0 3 
F8l4w 190.0 4 
2012-09-25 271.3 F606W 280.0 3 
F814W 200.0 4 
2013-05-13 99.9 F606W 280.0 3* 
F814W 200.0 4* 
2013-10-22 274.6 F606W 285.0 3 
F814W 285.0 4 
2014-10-26 275.2 F606W 265.0 3 
F814W 265.0 4 
2017-08-11 255.2 F606W 250.0 3* 
F814W 250.0 4* 
2017-08-29 268.3 F606W 250.0 3 
F814W 250.0 4 
2021-10-01 272.0 F606W 407.0 5 
F814W 307.0 6 





NotEe— Asterisk (*) denotes observations excluded from analysis. 


certainties, as well as the procedure used to determine 
them, are detailed in Appendix A. The additive error 
term dominates over the orpyg/VN term for bright 
stars, which makes the positional errors more uniform 
across epochs and magnitude as would be expected for 
systematic errors. 

As we are interested in the astrometry of the target, 
the target itself is not used to establish the transforma- 
tion into a common reference frame (§4.2.2) nor to judge 
the quality of the final transformation (Figure 9). 


4.2.2. Reference star selection 


Table 4. Gaia EDR3 Values 
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Parameter MB09260 MB10364 OB110037 OB110310 
Source ID 4064007633015639552 4042290560398692096 4056117808133831936 4056344036933003264 
RA (deg) 269.619073405 + 7.3e-07 269.272538687 + 1.9e-08 268.982636391+ 3.0e-08 267.855757518 + 2.3e-07 
Dec (deg) -26.839323825 + 5.9e-07 -34.451415987 + 1.5e-08 -30.561059419+ 2.5e-08  -30.409776355 + 1.7e-07 
Has (mas/yr) 7 -7.43 + 0.08 2.40 + 0.13 -2.08 + 1.12 
bs (mas/yr) = -6.80 + 0.05 -3.91 + 0.09 -6.75 + 0.58 
am (mas) = 0.40 + 0.08 0.15 + 0.13 0.54 + 1.16 
ZP-corrected 7 (mas) = 0.43 0.19 0.53* 
G (mag) 19.216 + 0.004 16.086 + 0.002 17.477 + 0.001 20.051 + 0.010 
RP (mag) = 14.929 + 0.009 16.323 + 0.010 - 
BP (mag) = 16.557 + 0.011 19.049 + 0.032 - 
ipd_gof_harmonic_amplitude 0.089 0.064 0.036 0.042 
ipd_frac_multi_peak 0 15 0 0 
ipd_frac_odd_win 18 0 0 55 
ruwe = 1.388 0.971 0.981 
astrometric_excess_noise (mas) 1.241 0.406 0.000 0.894 
astrometric_excess_noise_sig 2.657 12.020 0.000 0.332 
astrometric_params_solved 3 95 31 95 
phot_bp_rp_excess_factor = 1.69 1.39 za 





NotTe— Magnitude uncertainties are estimated from the Gaia reported flux errors. Zero-point (ZP) correction comes from Lindegren et al. 
(2021a), * denotes values that are extrapolations. For full descriptions we refer the reader to Gaia EDR3 documentation (van Leeuwen et al. 


2021), Section 13.1.1 gaia_source. OB110462 is not in Gaia. 


Reference stars are stars assumed to have linear proper 
motions, which are used to derive the reference frame 
transformation and the motions of the other stars. The 
selection of reference stars depends on multiple consider- 
ations, such as the stellar density, amount of geometric 
distortion, instrumental systematics, and number and 
brightness of targets of interest, to name a few. The goal 
is to balance having enough stars to establish the refer- 
ence frame, while excluding stars which would produce 
a non-stable reference frame. The criteria for reference 
star selection for each target is summarized in Table 5. 
We choose reference stars with brightness similar to the 
target, relatively large radial separations from our tar- 
get of interest, and exclude likely foreground stars. The 
target itself is also excluded from being a reference star. 
We detail the reasoning for these choices below. 

Brightness range: Due to the nature of CTE, there are 
strong magnitude-dependent astrometric residuals, even 
when using the latest CTE-corrected _flc images. How- 
ever, this is not unexpected (Kuhn & Anderson 2021). 
For this reason, stars are chosen to be in a brightness 
range similar to the one spanned by the target as nar- 
row as possible. For OB110310 and OB110462, which 
are relatively faint and where there are many stars of 
similar brightness, all stars falling within +0.1 mag of 
the target’s brightest and faintest in the HST data are 
used. MB09260, on the other hand, is brighter, with less 





stars of similar brightness, so the range is larger, with all 
stars falling within +1 mag of the target’s brightest and 
faintest in the HST data are used. MB10364 is so bright 
that many stars of comparable brightness are saturated 
in the longer exposures. Because of this, only bright 
(F814W < 18.0, F606W < 19.2) and unsaturated stars 
were selected; the bright limit on the magnitude range 
differs between epochs because of the different exposure 
times. 

Spatial separation: Only reference stars within 30”, 
or 20” for the denser field around OB110462, are used 
as reference stars. This minimizes the impact of geo- 
metric distortion residuals and spatially dependent PSF 
variations. 

Foreground stars: A key assumption in the astromet- 
ric alignment process is that reference stars have linear 
proper motions, and parallax effects are ignored. For a 
typical bulge star 8 kpc from Earth, this is a reasonable 
assumption as the parallax will be 1/8000 arcsec = 0.125 
mas, below our achievable astrometric precision. How- 
ever, for nearby stars, ignoring parallax is an issue when 
trying to derive an accurate transformation. As all the 
target fields are toward the highly extincted Galactic 
bulge, bright blue stars as identified on a CMD (Figure 
2) are likely to be nearby and have a non-negligible par- 
allax, and are excluded from the set of reference stars. 
The color-magnitude exclusion criteria are listed in the 
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last column of Table 5. For MB10364, no bright blue 
stars were removed as reference stars, as all the obser- 
vations came from within 6 weeks of the same time of 
year. Hence, any type of yearly parallax signal would 
be negligible within this time span. 

Number of detections: We require reference stars to 
be detected in most, if not all, epochs. If there are Nep 
total epochs, we require reference stars to be detected 
in Nep,detect = Nep — 2 epochs. 

Lastly, as the motion of the target is the quantity we 
are interested in, we do not use it as a reference star. 


4.2.3. Derived stellar proper motions 


To evaluate the goodness of the fits of the derived 
stellar proper motions, we consider the y? distributions 
of the position residuals 


an 
wD 


t 


where « = RA, y = Dec, (x,y); are the positions in 
the data, (x, y):, fie are the positions as derived from the 
linear motion fit, and o(z,y)4 are the positional uncer- 
tainties at time t. The distributions of x? and x7 for the 
reference stars detected in all epochs Nep are shown in 
Figure 9, with the expected x? distribution overplotted 
on top. The distributions for the positions in F814W 
and F606W are shown separately. The expected resid- 
ual distribution has Nep detect — 2 degrees of freedom, 
as there are two free parameters in the linear motion 
fit (initial position and proper motion). Note that un- 
like Gaia, this linear model fit does not include parallax. 
Parallax is only included when modeling the microlens- 
ing event. 


4.2.4. OB110462 bias correction 


There is a bright star (“the neighbor”) ~10 pixels 
(~0.4 arcsec) west of OB110462. The neighbor is ~3 
magnitudes brighter than OB110462 at baseline (F814W 
= 16.7 mag, F606W = 19.0 mag). Because of its proxim- 
ity and high contrast, the neighbor’s PSF might “leak” 
onto OB110462 and alter its astrometry and photome- 
try. We perform injection and recovery tests to ascertain 
the reliability of faint source extraction near a bright 
source, in order to determine whether the astrometry 
and photometry of OB110462 as determined in §4.2 and 
84.3 is biased by the bright neighbor. 

The methodology and results of the injection and re- 
covery analysis are detailed in Appendix B. In summary, 
the positional bias is negligible when the target is highly 
magnified and of similar brightness to the neighbor. 


However, in epochs where the target is no longer mag- 
nified, the bright star biases the position of the target. 
In F814W, where the resolution is lower, the measured 
position of the target is biased toward the neighbor by 
~ 0.4 mas along the target-neighbor separation vector. 
In F606W, where the resolution is higher, the bias is less 
(~ 0.25 mas) with the direction of bias more randomly 
oriented. Similarly, the photometric bias is larger when 
the contrast is large, with the bright neighbor causing 
the extracted photometry of OB110462 to be brighter 
than the injected values. The effect is again more severe 
in F814W than F606W because of the lower resolution. 

Using the results of the injection and recovery analy- 
sis, we calculate a bias correction to apply to OB110462 
astrometry and photometry (Table 17). The values in 
the table are added to astrometry and photometry de- 
rived in §4.2 and 84.3; the uncertainties are added in 
quadrature to the uncertainties in §4.2 and §4.3. 

We only perform this analysis for OB110462, as it is 
the only faint target near a bright companion. All the 
other targets are either bright with faint companions, 
isolated, or both bright and isolated. 


4.2.5. Astrometric color offset 


As mentioned in §4.1 the data taken in F606W and 
F814W filters are treated as independent measurements. 
For OB110037 and OB110462, the astrometric measure- 
ments in F606W and F814W do not agree within the 
uncertainties. 

For OB110037, although the 2011 and 2012 epochs 
show good agreement, the 2013 to 2017 epochs become 
increasingly discrepant as time goes on. We attribute 
this difference to binarity (87.4). 

In contrast, for OB110462 the astrometry in the 
F606W and F814W are discrepant in all datasets, but 
the difference appears to be a relatively constant offset 
with time. This is true both before and after applying 
the bias correction in §4.2.4. Because the nature of the 
color difference appears to be a constant offset, we apply 
a constant shift to the OB110462 F606W astrometry 








_ ‘3 Wa,t(LFEsiaw,t _ LFGO6W,t) 
Ke (9) 
ye Wet 
oy Wyt (Yrsi4w,t — YFooew,t) 
Ay = (10) 
Die Wy,t 
where 
Wet = (02 rsiaw.t + o2 reoswt) (11) 
Wyt = (07 peiaw.t + a ce) (12) 
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Figure 9. Histogram of x” residual values to the linear fits with no parallax for the reference stars of each target (Table 5). 
In each panel, the left column shows the distributions for reference stars in F606W, while the right column shows F814W. The 
number of reference stars is listed as N. In each panel, the top row shows the ? distribution of residuals of a linear fit to 
positions vs. time in x = RA and y = Dec. The bottom row shows the y? distribution of residuals of a constant fit to magnitude 
vs. time. The expected x? distributions are shown in black, with the number of degrees of freedom listed as Nao Pe 


and t indexes the observation times. Thus the modified 
astrometry for F606W is 


(13) 
(14) 


y 
Treoew,t = TFeoew,t + Ax 


/ 
Yreosw,t = YFeoew,t + Ay. 


The values of the offset are A, = —0.57 mas and A, = 
0.39 mas. Note that these offsets are calculated using the 
bias-corrected astrometry.? See Appendix H for further 
details and justification. 

We also investigate stars nearby to determine whether 
any of them show similar behavior. For the 70 stars 
within 3” of OB110462, we calculate the average posi- 
tional offset between F814W and F606W in RA and Dec 
using Eqs. 9 and 10. We then search for stars where the 
average positional uncertainty in F814W and F606W 
(whichever is larger) is smaller than the average posi- 


9 Note that even before the bias correction of §4.2.4 is applied, 
this color offset is still present. In fact, it is slightly larger, with 
Az = —0.79 mas and A, = 0.52 mas. The bias correction is not 
the source of the color dependent astrometric offset; rather, it 


helps to slightly decrease the offset. 


tional offset to determine which differences are signif- 
icant. There are 4 stars where the average positional 
offset is greater than the average positional uncertainty 
in RA, and an additional 4 stars the average offset is 
greater than the average positional uncertainty in Dec. 
Hence, a total of 8 out of the 70 stars near the target 
also show these significant offsets. Thus, this effect is 
seen for roughly 10% of the stars, and so is not very 
unusual. Although we currently have no explanation for 
its significance, it appears random, and thus include the 
astrometric offset when analyzing the data, although it 
may be attributed to binarity (see §7.4 and §4.1 in the 
main paper). See Appendix H for additional details. 


4.2.6. Comparison to Gaia proper motions and parallazes 


For the three targets with astrometric solutions in 
Gaia EDR3 (MB10364, OB110037, and OB110310; Ta- 
ble 4), the Gaia proper motions and parallaxes are com- 
pared to the fit proper motions and parallaxes presented 
in §7 and Tables 11 - 13. Direct comparisons are made in 
Table 6 and Figures 10 and 11. Note that we fit a proper 
motion and parallax to the source and lens along with 
an astrometric microlensing model. On the other hand, 
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Table 5. 


Reference star criteria 








Target 


Magnitude range 


Radius  Nep, detect 


EXCLUDED bright blue stars 





MB09260 





Target range + 1 mag: 
14.5 < F814W < 18.8 
16.9 < F606W < 21.7 


30” 11 


F606W-F814W < 2.6 
& F814W < 20.6 





MB10364 


2010-09-13: 

12.5 < F814W < 18.0 
11.5 < F606W < 19.2 
2010-10-26: 

15.1 < F814W < 18.0 
13.4 < F606W < 19.2 
2011-10-31: 

15.8 < F814W < 18.0 
14.4 < F606W < 19.2 
2012-09-25: 

16.0 < F814W < 18.0 
14.5 < F606W < 19.2 
2013-10-24: 

16.1 < F814W < 18.0 
14.7 < F606W < 19.2 


30” 8 


None excluded 





OB110037 


Target range + 0.5 mag: 


14.4 < F814W < 16.9 
16.4 < F814W < 18.8 


30” 12 


F606W-F814W < 1.75 
& F814W < 19.6 





OB110310 


Target range + 0.1 mag: 


16.9 < F814W < 18.7 
19.6 < F606W < 21.4 


30” 12 


F606W-F814W < 2.4 
& F814W < 21.0 








OB110462 





Target range + 0.1 mag: 





17.1 < F814W < 20.0 
19.2 < F606W < 22.1 


20” 14 


F606W-F814W < 1.9 
& F814W < 20.6 








NOTE 


These are the criteria for the last pass. 


Table 6. Gaia vs. HST Proper Motions 



































Target HusT,L (mas/yr) LuST,s (mas/yr) UGaia (mas/yr) 
MB10364 = (—5.11 74-62, —7.787 0-98) (—7.5670:12, —6.497?-11) (—7.43 + 0.08, —6.80 + 0.05) 
OB110037 (6.275139, —6.56To:1) (2.195534, 3.877929) (2.40 + 0.13, —3.91 + 0.09) 
OB110310 (—0.021}-93, —4.6813°79) (—2.417915, —7.2619°08) (—2.08 + 1.12, —6.75 + 0.58) 














NotTe—The source and lens proper motions here have been transformed into the absolute Gaia proper 
motion frame, which is offset to the HST proper motion frame as described in Appendix C. The 
uncertainties on wy, and jg also reflect the uncertainty in the Gaia to HST proper motion trans- 
formation; the standard error on the mean of that transformation was added in quadrature to the 
uncertainties from the proper motion fits. For this reason, the uncertainties for uz and ug in this 


table do not match those in Tables 11 - 13. 
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Figure 10. Proper motions for stars in the field of the target. The proper motions derived from F814W HST observations for 
stars within 30” of the target are shown in gray with lo uncertainties. Only stars with F814W < 23 for MB09260, OB110310, 
OB110462 and F814W < 22 for MB10364, OB110037 are shown. For OB110462 there are two models depending on the 
likelihood used (default weighted “DW” or equal weighted “EW”, see §5.3 for details). As the alignment procedure places the 
stars in a reference frame where the relative motion is zero, a constant offset must be added to obtain proper motions in the 
original Gaia absolute reference frame. This offset is calculated by matching the stellar positions in HST to those in Gaia with 
astrometric_excess_noise_sig< 2, and then calculating the 3c clipped average weighted by the uncertainty in their difference. 
The offset value is given in the title for each field. For targets in Gaia where a single-star proper motion is estimated (MB10364, 
OB11037, OB110310), they are plotted as black squares. The lens and source proper motion as determined from fitting the 


HST data with a microlensing model is shown in blue and red 1 — 2 — 30 contours, respectively. Note that the red contours are 
extremely small. 
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Figure 11. Comparison of Gaia parallaxes with the fit lens 
and source parallax. 


Gaia fits a proper motion and parallax to the source and 
lens (the “target” ) as they are unresolved, and assumes 
the target is a single star with parallax. The effect of 
astrometric lensing in Gaia on the proper motions is 
negligible (Appendix D). Additional details about the 
various Gaia metrics discussed are in Appendix E and 
Lindegren et al. (2021a). 

To make proper motion comparisons between Gaia 
and HST, the proper motions from the HST frame in 
which the fitting was performed need to be transformed 
into the absolute Gaia frame, as the iterative astromet- 
ric alignment procedure described in §4.2 produces a 
reference frame that is at rest with the average proper 
motion of the aligned stars (Lu et al. 2016). See Ap- 
pendix C on how the proper motion offset between the 
HST and Gaia frames is calculated. The source and lens 
proper motion in the Gaia frame, as well as the Gaia tar- 
get proper motions, are listed in Table 6 for MB10364, 
OB110037, and OB110310. The vector point diagram 
for all stars from the HST observations transformed to 
the Gaia frame, along with the 1 — 2—30 contours from 
the source and lens fit, are shown in Figure 10. For 
MB10364, OB110037, and OB110310 the target proper 
motion from Gaia is also included. 


MB10364—MB10364’s proper motion in Gaia is 
(—7.43 + 0.08, —6.80 + 0.05) mas/yr. The fit from HST 
for the lens is (—5.1143:94, —7.7843-38) mas/yr, incon- 
sistent with Gaia in RA and Dec at ~ 2c. The fit from 
HST for the source is (—7.567$'13, —6.4972: 51) mas/yr, 
inconsistent with Gaia in Dec at ~ 20. MB10364 has 
parallax 7 = 0.43 + 0.08 in Gaia. The source and lens 
parallax from the MB10364 fits are mg = 0.111003, 
T= O24 TR. neither of which are consistent with the 
Gaia value. The Gaia fit for MB10364 has a large renor- 
malized unit weight error (RUWE = 1.388), and a large 
astrometric excess noise (€ = 0.406 mas) with high sig- 
nificance (D = 12.020), indicating the single star model 
is not providing a good fit. This mismatch is not due to 
astrometric microlensing (Appendix D). The most likely 
explanation for the discrepancies are crowding-— there are 
several stars close to MB10364 that would cause confu- 











sion (Figures 1 and 12). As the source and lens are not 
resolvable, the fact that 15% of the Image Parameters 
Determination (IPD) algorithm has identified a double 
peak is likely due to confusion. This means the Gaia 
measurement is not reliable. See Appendix E for fur- 
ther discussion. 


OB110037—OB110037’s proper motion in Gaia is 
(2.4040.13, —3.91£0.09) mas/yr. The fit from HST for 
the lens is (6.277 }'39, —6.56104,) mas/yr, inconsistent 
with Gaia in both RA and Dec at ~ 30. The fit from 
HST for the source is (2.197334, —3.8710-20) mas/yr, 
consistent with Gaia. The caveat to this is that the 
fit to the HST F606W astrometry is poor (Figure 18). 
However, the astrometric lensing model is mostly inde- 
pendent of the proper motion model, hence the proper 
motion value can still be believed. OB110037 is very 
well measured and behaved in Gaia (it is the only target 
with a 5-parameter solution, see Table 4). Additionally 
from the image (Figures 1 and 12) OB110037 is rela- 
tively bright and isolated. The fact that the fit source 
and Gaia proper motions are consistent would indicate 
that the lens is dim in comparison to the star in Gaia 
G band. This makes sense as the astrometric shift for 
OB110037 is undetectable by the time Gaia begins ob- 
serving. 

OB110037 is well modeled in Gaia by a single source 
with parallax. Although the microlensing model is a 
poor fit to the astrometry due to a time-dependent 
color-offset, the proper motions from our model are 
in good agreement with Gaia. OB110037 has paral- 
lax + = 0.19 + 0.13 in Gaia. The source and lens 
parallax from the OB110037 fits are ms = 0.121305, 
my = 0.587913. The source parallax is consistent with 
the Gaia value. We note that the source may appear 
well-behaved in Gaia because the astrometry is in a sin- 
gle filter. The multi-band HST astrometry may be use- 
ful in identifying binary companions. See Appendix E 
for further discussion. 











OB110310—OB110310’s proper motion in Gaia is 
(—2.08 + 1.12, —6.75 + 0.58) mas/yr. The fit from 
HST for the lens is (—0.02*}-92, —4.6813-73) mas/yr, 
consistent with Gaia due to large uncertainties in both 
Gaia and the fit. The fit from HST for the source is 
(—2.4172:53, —7.2610-08) mas/yr, consistent with Gaia. 
OB110310 has parallax + = 0.53 + 1.16 in Gaia (note 
the OB110310 zero point correction is an extrapola- 
tion), consistent with non-detection of parallax. The 
source and lens parallax from the OB110310 fit are 
mg = 0.107303, m, = 0.221957, which are both con- 
sistent within the very wide uncertainties of Gaia. The 
fact that the source and Gaia proper motions are consis- 











MB09260 
T T 


MB10364 
T T 


0B110037 
T T 


17 


0B110310 0B110462 
T T T T 





Dsrr 





0.2- : 44 4h 

















0.0- 1 ! n Io ! 1 ! JIC 





DsFF,max,Fe14w | 
* Bser, max, F606W 
Dsre, rit J 
Typical seeing 


























t 
2.0 0.5 


\ 
1.0 


1 
1.5 2.0 


Seeing disk radius (arcsec) 


Figure 12. Comparison of bsrr inferred from fitting microlensing events, as compared to upper limits inferred from the 
high-resolution HST images. The typical seeing disk radius ranges are shown in light gray for MOA (MB09260, MB10364) and 
OGLE (OB110037, OB110310, OB110462). The inferred value of bsrr from fitting each microlensing event is shown as a black 
line, with uncertainties in dark gray. The upper limits on blending bsrrmaz as a function of aperture radius, as estimated 
using the method described in §7.7, are shown as the solid red line (F8&14/W) and the dotted blue line (F606W). For the three 
OGLE targets, at the relevant seeing disk radii, the inferred value of bsrr from the fit are at or below bsrrmaz as inferred 
from the HST images. For the two MOA targets, the inferred values of bsrr are higher than the maximum value as inferred 
from the HST images. In the case of MB09260, this is likely due to the bright star around 1.0” from the target which causes 
an abrupt change in bsrr around that radius; proper convolution with a PSF would likely resolve the difference. In the case of 
MB10364, the difference is minimal, and likely just due to the imperfect nature of the comparison (for example, lack of proper 
PSF convolution, and F814W not exactly matching the red filter in MOA/OGLE). 


tent would indicate that the lens is dim in comparison 
to the star in Gaia G band; blending in F814W and 
F606W also suggest a dark lens (Table 13). The fact 
that in the IPD 55% of transits have either truncation 
or multiple gates flagged in one or more windows indi- 
cates likely contamination. OB110310 is not very bright 
and in a somewhat crowded region (Figures 1 and 12). 
The astrometric noise is large (« = 0.894 mas) but the 
value is insignificant (D = 0.332). Together these ex- 
plain why the Gaia measurement does not produce very 
good constraints. See Appendix E for further discussion. 


4.3. HST Photometric Analysis 


To obtain precise relative photometry, for each epoch 
a small constant magnitude offset is applied to the stars. 
The offset is calculated by assuming the reference stars 
have constant brightness in time, which we define as the 
30 clipped mean. As with the positional uncertainties, 
for the magnitude uncertainties we used the error on 
the mean, with an additional additive error empirically 
determined during the astrometric alignment process; 
details are in Appendix A. 

Analogous to the positional transformation, we evalu- 
ate the magnitude transformation by checking how well 
a constant magnitude describes the stars 


oi : 


where m, are the calibrated magnitudes, mg is the con- 
stant magnitude fit, and om, is the uncertainty on the 


ile a) (15) 


Om; 


calibrated magnitude. The x? distributions for the mag- 
nitude residuals of the reference stars detected in all Nep 
epochs are shown in the bottom panels of Figure 9. The 
expected residual distribution has Ne» detect — 1 degrees 
of freedom, as there is one free parameter in the constant 
magnitude fit. 


5. MICROLENSING MODELING 


The mass of the lens M; in a microlensing event is 
given by 


(16) 


where 0 is the angular Einstein radius (Equation 17), 
wm is the microlensing parallax (Equation 21), and 
4G/(1AU - c?) = 8.14mas/Mo is a constant. 
Densely sampled photometric microlensing observations 
constrain tz, while astrometric microlensing observa- 
tions constrain 6,.1° 

To measure these quantities of interest, we simulta- 
neously fit the ground-based photometry and HST pho- 
tometry and astrometry with a point-source point-lens 
(PSPL) microlensing model with parallax. We do not 
consider models involving either binary lenses or sources, 
nor higher-order effects beyond parallax; these are be- 
yond the scope of this work. Discussion of the need 


K 


10 Theoretically, astrometric microlensing observations should also 
be able to constrain 7g, but due to the cadence of observations, 


this is currently unachievable. 
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Table 7. HST Calibrated Data for Each Target 

















































































































































































































































Target Filter Date RA (mas) Dec (mas) Mag (Vega) 
MBO09260 F814W 2009-10-19 7.80 + 0.15 5.15 + 0.14 15.484 + 0.005 
2010-03-22 5.31 3.80 17.656 + 0.005 
2010-06-14 4.31 2.97 17.760 + 0.007 
2010-10-20 2.544 1.60 4 17.812 + 0.014 
2011-04-19 0.12 + 0.01 + 17.838 + 0.018 
2011-10-24 -2.44 -1.71 17.833 + 0.005 
2012-09-25 -7.07 -4.80 17.836 + 0.007 
2013-06-17 -10.57 40.26 -7.02 + 17.829 + 0.011 
F606W 2009-10-01 8.15 5.57 4 17.899 + 0.013 
2011-04-19 0.13 0.17 20.750 + 0.083 
2011-10-24 -2.17 -1.59 20.748 + 0.007 
2012-09-25 -6.23 + -4.41 4 20.733 + 0.023 
2013-06-17 -10.10 + -7.21 + 20.738 + 0.036 
MB10364 F814W 2010-09-13 9.61 8.24 13.366 + 0.011 
2010-10-26 14.657 + 0.005 
2011-10-31 15.315 + 0.006 
2012-09-25 15.315 + 0.010 
2013-10-24 15.316 + 0.009 
F606W 2010-09-13 14.538 + 0.017 
2010-10-26 15.842 + 0.017 
2011-10-31 16.498 + 0.006 
2012-09-25 16.504 + 0.008 
2013-10-24 16.502 + 0.008 
OB110037  F814W ~— 2011-08-15 14.864 + 0.021 
2011-09-26 15.029 + 0.005 
2011-11-01 15.774 + 0.006 
2012-05-07 16.315 + 0.019 
2012-09-25 16.321 + 0.009 
2013-10-21 16.327 + 0.005 
2014-10-26 16.328 + 0.006 
2017-03-13 16.334 + 0.007 
2017-09-04 16.322 + 0.011 
F606W 2011-08-15 16.916 + 0.011 
2011-09-26 17.086 + 0.020 
2011-11-01 17.794 + 0.005 
2012-05-07 18.301 + 0.016 
2012-09-25 é 18.306 + 0.006 
2013-10-21 0. 18.314 + 0.006 
2014-10-26 3. 18.320 + 0.015 
2017-03-13 9. 18.331 + 0.008 
2017-09-04 9. 18.326 + 0.027 
OB110310 F8l4W 2011-09-21 5. 16.945 + 0.013 
2011-10-31 4. 18.058 + 0.004 
2012-04-24 3. 18.602 + 0.013 
2012-09-24 2. 18.616 + 0.011 
2013-10-21 0. 18.621 + 0.005 
2017-03-14 _-7.92 18.608 + 0.005 
2017-09-01 -8.88 18.613 + 0.012 
F606W 2011-09-21 5.43 19.663 + 0.018 
2011-10-31 4. 20.780 + 0.007 
2012-04-24 3. 21.329 + 0.034 
2012-09-24 2. 21.180 + 0.212 
2013-10-21 21.339 + 0.007 
2017-03-14 -7.81 4 21.333 + 0.022 
2017-09-01 _-9.31 21.335 + 0.034 
OB110462.  F814W 2011-08-08 7.53 17.209 + 0.028 
2011-10-31 6.44 4 18.849 + 0.006 
2012-09-09 4.08 + 19.756 + 0.009 
2012-09-25 4.25 19.767 + 0.008 
2013-10-22 1.43 19.839 + 0.048 
2014-10-26 = -0.87 + 19.881 + 0.009 
2017-08-29 -7.18 + 19.874 + 0.009 
2021-10-01 -15.67 -24.93 19.865 + 0.020 
F606W 2011-08-08 6.76 11.77 19.313 + 0.023 
2011-10-31 6.25 + 10.54 + 20.974 + 0.010 
2012-09-09 4.60 + 5.87 + 21.867 + 0.050 
2012-09-25 4.38 6.40 21.920 + 0.010 
2013-10-22 1.76 2.45 22.011 + 0.015 
2014-10-26 = -1.22 + -0.83 4 22.042 + 0.015 
2017-08-29 -7.10 + -10.62 + 22.017 + 0.027 
2021-10-01 -15.44 -25.06 22.021 























Norr— Relative positions and magnitude of the target by epoch and filter. 


for models more complex than PSPL with parallax are 
discussed in 87. 

Throughout this section and the remainder of the pa- 
per, we report vector quantities decomposed into their 
RA and Dec components, subscripted by “ax” and “6”. 
Because we work solely in the Equatorial coordinate sys- 
tem, we also equivalently refer to RA as East and Dec as 
North, subscripted by “E” or “N”, where RA increases 
to the East and Dec increases to the North. 


5.1. Microlensing Definitions 


All microlensing quantities defined in the following 
section are in the heliocentric reference frame. 

By rearranging the terms in Equation 16, the Einstein 
radius, which sets the angular scale of the microlensing 
event, can be written 


On =v KT re Mrz (17) 


where 7,¢] = Ty — Tg is the relative parallax of the lens 
and source. 

The Einstein crossing time tg, the time it takes for 
the source to traverse the angular radius of the lens and 
sets the timescale of the events, is given by 

Op 


i= (18) 
Lrel 





where /t;re, is the lens-source proper motion |pg — p17]. 

The source-lens separation on sky 05 — 0; normalized 
by the Einstein radius is denoted u(t). The minimum 
separation is denoted u(to) = uo. The impact param- 
eter uo is the scalar quantity associated with uo. We 
follow the convention of Gould (2004) where if uo,¢ > 0, 
the source is to the East of the lens, and up > 0; if 
uo,E < 0, the source is to the West of the lens, and 
ug < 0.1! If rectilinear motion is assumed, the lens- 
source separation is given by 


u(t) =u + Ayer (19) 
E 

However, an Earthly observer’s perspective of the lens- 

ing event is modulated by the Earth’s motion around the 

Sun. For events with long duration (tg 2 3 months), the 

Earth’s orbital motion violates this rectilinear assump- 

tion and must be taken into account. This modifies u(t) 


11 Note that we define our coordinate system differently than Gould 
(2004). Gould (2004) works in a system where the position of the 
source relative to the lens is defined in a coordinate system that is 
right-handed in the relative proper motion and minimum separa- 
tion vector. However, we work in a system where the coordinate 
system is consistent on the sky; this means a coordinate system 
based on the relative proper motion and minimum separation 


vector does not always preserve handedness. 
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to 
t—to. 
u(t) = Uo + SF Piret ~~ TpP(t) (20) 
where ae 
m= (21) 


is the microlensing parallax and P(t) is the parallax 
vector, defined to be the Sun-Earth separation vector 
normalized by 1 AU. The microlensing parallax vector 
aw (not to be confused with the parallax vector P(t)) 
encodes the magnitude of the microlensing parallax and 
the direction of the relative source-lens proper motion: 


TE = TE Pret: (22) 
The photometric brightening of the source is given by 


Ai —. (23) 


where the total flux F(t) in the telescope aperture is 
F(t) = A()Fs + Fi + Fw (24) 


where F's, F’,, and F'y are the fluxes of the source, lens, 
and neighboring un-lensed stars in the aperture of the 
telescope, respectively. The source flux fraction is 
Fs 

a Fs + Fy + Fn ci 
and quantifies the fraction of light lensed in an observed 
microlensing event. Note that bgrr depends on the 
observing wavelength and seeing/aperture. The non- 
source flux Fy, + Fy is also called blend flux. Blend flux 
decreases bsp and dilutes the magnitude of both the 
photometric amplification and astrometric shift. 

A PSPL photometric microlensing event is character- 
ized by five geometric parameters: to, uo, tg, and mz. 
For each telescope that observes this event, two addi- 
tional parameters, Mbase and bgrr, are needed to de- 
scribe each lightcurve; these depend on the seeing and 
camera filter. As all photometric microlensing observ- 
able quantities are normalized by @z, the lens mass can- 
not be determined. 

In contrast, astrometric microlensing provides a direct 
measurement of 07. The apparent position of the source, 
i.e. the centroid of the lensed source images 05.., is given 
by : 

(u* + 3)uOe 

95 c(u, Oz) = ro aa 

Assuming no blended light, the difference between the 

source’s apparent and true positions, i.e. the astrometric 
shift, is given by 


(26) 


ub rp 


6-(u, Oe) = wy’ 


(27) 
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The astrometric shift 6. is maximized at u = /2. This 
corresponds to the value of the maximum astrometric 
shift 
ee 
cman — V3" 


which is directly proportional to the Einstein radius. 

Note that Equation 27 is for an unblended event, i.e. 
bsrr = 1. If an event has bgsrr < 1, the non-source 
light would dilute the astrometric shift (see Equations 
11 and 12 in Lam et al. (2020)). For the five candidates 
analyzed here, bgrr ~ 1 in the HST filters, so assuming 
that the astrometry is unblended in HST is valid. 


(28) 


5.2. Modeling framework 


We perform parameter estimation using a Bayesian 
framework. Bayes’ theorem 


™(O)L(y|O) = Z(y)P(Oly) (29) 


relates the prior 7 and likelihood £ to the evidence Z 
and posterior P. The goal of parameter estimation is to 
calculate P. The likelihood L(y|@) is presented in Ap- 
pendix F, and the priors 7(@) are discussed in Appendix 
G. 

The data is fit using MultiNest (Feroz et al. 2009), 
an implementation of the nested sampling algorithm 
(Skilling 2004). Nested sampling produces an estimate 
of the evidence Z = f[ 7(@)L(O)dO, and as a by- 
product, the posterior P. In contrast to methods such 
as MCMC, nested sampling is designed to better explore 
multimodal likelihood spaces; however, care must still be 
taken to ensure that all local minima are explored. 


5.3. Likelihood weighting 


There is a question of how best to combine the pho- 
tometry and astrometry data sets, as they represent two 
different types of measurements. In particular, the ques- 
tion is how much weight each dataset should receive 
given that there are several orders of magnitude more 
ground-based data points than HST data points. Thus, 
the ground-based photometry has an outsize effect on 
the likelihood. It can be argued that each data point 
should contribute equally to the likelihood. We con- 
sider this to be “default weight” (hereafter abbreviated 
as DW) likelihood, i.e. 


log Lot = log LO, phot =f log L phot =F log Litast (30) 


where Lo phot, £H,phot, and Ly,as¢ are the likelihoods of 
the OGLE or MOA photometry, HST photometry, and 
HST astrometry respectively. 

However, it could also be argued that each dataset is 
independent, and so should each contribute equally to 


the likelihood. We call this the “equal weight” (hereafter 
abbreviated as EW) likelihood, i.e. 


log LO,phot 4 log Lit phot Aly log Linast ( 
NO, phot 


log Ltot = 31) 








NH, phot NH, ast 


where NO,phot, NH,phot; and Ny,ast are the number of 
data points in in the OGLE photometry, HST photom- 
etry, and HST astrometry. 

These different likelihoods are essentially giving the 
different datasets different weights, in the case that they 
are inconsistent with each other. For all targets, we fit 
using the default weight likelihood; for OB110462 we 
additionally fit using the equal weight likelihood. 


6. CONSTRAINING THE NATURE OF THE LENS 


By modeling photometric and astrometric microlens- 
ing data as described in §5, the lens’ mass can be mea- 
sured. However, a mass measurement alone cannot dis- 
tinguish the difference between a 5Mo stellar lens and a 
5Mo BH. Additional information about the lens’ bright- 
ness is needed, which can be estimated using the source 
flux fractions in the high-resolution HST filters. By cal- 
culating the brightest star allowed by the inferred source 
flux fractions, we can determine whether a luminous lens 
(i.e. stellar lens) or dark lens (i.e. compact object lens) 
scenario is more likely. 

We follow a procedure similar to Wyrzykowski et al. 
(2016) in order to calculate the probability of a dark 
lens. The two main differences between the analysis of 
Wyrzykowski et al. (2016) and ours is that they work 
with photometric data only and must invoke a Galactic 
model to obtain lens masses and distances, while in our 
case astrometric data allows us to fit lens mass and dis- 
tances. This greatly simplifies the analysis as we do not 
need to calculate the Jacobian to transform the PDFs 
between different variables. In addition, we use realistic 
stellar and Galactic models to determine the lens lu- 
minosities, instead of relying on simple mass-luminosity 
scaling relations. 

From fitting the microlensing data, we have posterior 
distributions for the lens mass Mz, distance dz, baseline 
magnitude Mpbase, and source flux fraction bsrr. In the 
following analysis we reassign any fit values with bsrr > 
1 to instead be equal to 0.999999. This is because values 


of bsrr > 1 would result in an infinite lens magnitude 
(Equation 32).1? 

We draw a random sample of 1000 stars from our pos- 
terior. For each star we calculate 1) the brightest lens 
allowed by bsrr and Mpase, and 2) the brightest star 
allowed by My and dz. 

To calculate 1), by assuming there are no contami- 
nating neighbor stars, bspr = Fs/(Fs + F,), an upper 
limit can be placed on the brightness of the lens: 


ML = Mpbase — 2.5 logy 9(1 _ bsrr). (32) 


We denote this ML (bsrF, Mbase): 

To calculate 2), we use the simple stellar popula- 
tion synthesis code SPISEA (Hosek et al. 2020) to gen- 
erate a suite of simple stellar populations (SSPs) to 
simulate the possible lens population. We use the 
MISTv1.2 solar metallicity isochrones (Choi et al. 2016), 
get_merged atmosphere atmosphere model!*, Damineli 
et al. (2016) reddening law, and Kroupa (2001) initial 
mass function (IMF) over the mass range 0.1Mo < 
M < 120Mo. Each cluster is 10+*Mo, in order to re- 
duce stochastic effects in the sampling of the IMF. 

SSPs are generated at the distances spanned by each 
target’s lens distance posteriors, sampled every 0.25 kpc. 
SSPs of ages 7.0 to 10.0 log,g years in increments of 0.5 
log;, years are simulated at each sampled distance in 
order to cover the age range of disk and bulge stars. 
The stellar age distribution as a function of distance is 
drawn from the Galaxia Milky Way stellar simulation 
(Sharma et al. 2011), which implements a version of the 
Besancon Galactic model of Robin et al. (2003). As 
our target sight lines are toward the bulge, most of the 
stars are old, with = 85% of stars being at least 9.0 
log;, years old, although younger stars tend to be at 
closer distances as disk stars dominate. To calculate the 
simulated stars’ apparent magnitudes, we use the 3-D 
extinction map of Marshall et al. (2006), accessed via 
the dustmaps software package (Green 2018), to obtain 
an Ax, value for each distance. 

For each sample from the posterior, we sample one 
of the SSPs at the corresponding distance at a particu- 
lar stellar age from the Galactic distribution. All stars 


12 One might worry this could bias the results. However, a visual 
check comparing the posteriors (both individual and joint) for 
lens mass and distance don’t show a change when excluding fits 
with bsrr > 1. In particular, the lens mass and distances are 
nearly independent of the lens magnitude. In addition, if the 
true value of bgrr is 1, we would expect some scatter around 
that value, including unphysical > 1 values. 

13 This is a combination of the ATLAS9 (Castelli & Kurucz 2004), 
PHOENIX v16 (Husser et al. 2013), BTSettl (Baraffe et al. 2015), 
abd Koester10 (Koester 2010) models. For further details see 
Appendix B of Hosek et al. (2020). 
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in the simulated SSP with masses within 10% of the 
sampled posterior mass are identified, and the brightest 
apparent magnitude of the star in this group is denoted 
mz(M,, Dr, *). If no stars are found within 10% of 
the sampled posterior mass, this indicates that stars of 
that mass have all died and formed compact objects. 

If mri(Mz, Dr, *) a mr(bsrr,Mbase); a star with 
the inferred stellar mass from stellar evolution models 
would be too bright to be hidden in the blended light 
allowed by the fit. This means that the lens is not a star. 
We take this to mean the lens is a compact object, and 
hence a dark lens. Thus the lower limit on the probabil- 
ity of a dark lens in the observed filter is the fraction of 
samples where mzp(Mz,, Dri, *) < mz(bsrr,Mbase) OF 
no mass match is found. 

Objects with mz(M,, Dr, *) > mi (bsrr,Mbase) 
are samples where a luminous stellar lens is consis- 
tent with the inferred amount of blending. Note that 
mr(bsrr,Mbase) is an upper limit on the brightness of 
the lens, as we assumed all the blended flux in Equa- 
tion 32 is due to the lens; it could be due to unre- 
solved unrelated neighbor stars. Thus the upper limit 
on the probability of a stellar lens in the observed filter 
is thus the fraction of samples where mp(Mz,, Di, *) > 
mr (bsFr,Mbase): 

For dark lenses (i.e. any objects where 
mr(Mrz, Dr, *) < mr (bsrr,Mbase) or no star of the 
same mass was found in the simulated SSP), we catego- 
rize them as brown dwarfs (BD), white dwarfs (WD), 
neutron stars (NS), or BHs by simplistically sorting 
them by their masses: 


e Brown dwarfs: M < 0.2Mo 

e White dwarfs: 0.2Mo <M <1.2Mo 
e Neutron stars: 1.2Mo < M < 2.2Mo 
e Black holes: M > 2.2Mo. 


In reality there is overlap between white dwarf and neu- 
tron star masses, and the overlap between neutron star 
and BH masses (if they overlap at all) is unknown. In 
addition, the maximum brown dwarf mass set by stellar 
physics is 0.08M,, but we extend this up to 0.2Me to 
have continuity between the lowest mass WDs of around 
0.2Ma. Hence these values are only approximate. 

The above analysis is performed for both the HST 
F606W and F814W filters. The reported probabilities 
are the joint constraint. A lens is dark if no masses 
consistent with a stellar lens are found in either filter. 
A lens is also dark if the maximum inferred lens flux 
is insufficient to hide a star (ie. mz(Mz,Dzi,*) < 
mr(bsrr,Mbase)) in either filter. A lens is luminous if 
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Table 8. 99.73% credible intervals/upper limits 




















Parameter Om (mas) Se,maa (mas) Mp(Mo) a 
MB09260 < 2.42 < 0.85 a7? 72 o.09F 0-18 
MB10364 < 1.76 < 0.62 o.21F0-8% o.27t 0-02 
OB110037 1.247630 0.4470-39 0.4179-37 — o.3740-02 

$2.98 +0:20 
OB110310 < 2.75 < 0.97 o.7st2-88 0.13838 
OB110462 (BW) 4.13 72°05 1.461405 2.15t3-59  g.2at0-t8 
©B110462 (DW) 3.897369 1371000 3.791244 0,127 0-08 





Nore—There are two entries for OB110462, one based on a fit with “default 
weighting” (OB110462 DW) and one based on a fit with “equal weighting” 
(OB110462 EW). Sce §5.3 for details. 


Table 9. Lens type probabilities (%) 











Target Star BD WD NS BH 
MB09260 4 0 38 44 14 
MB10364 36 29 35 0 i) 
OB110037 74 0 26 0 i) 
OB110310 5 3 65 22 5 
OB110462 DW i) 0 0 0 100 
OB110462 EW 0 0 6 50 44 








Nore— The Star probabilities are upper limits, while 
the brown dwarf (BD), white dwarf (WD), neutron 
star (NS), and BH probabilities are lower limits. 
This is because the luminous lens probabilities are 
upper limits, and the dark lens probabilities are 
lower limits; see §6 for details. Note that there are 
two entries for OB110462, one based on a fit with 
“default weighting” (OB110462 DW) and one based 
on a fit with “equal weighting” (OB110462 EW). 
See §5.3 for details. 


the maximum inferred lens flux is sufficient to hide a 
star (i.e. mr(Mz, Dr, *) > mr (bsrr,Mbase)) in both 
filters. We do not perform this analysis for the OGLE 
or MOA photometry parameters as the high resolution 
HST images show the seeing-limited apertures have un- 
related neighbor stars in the blend, and hence the limits 
will all be weaker than using HST. 


7. RESULTS IN DETAIL 


A large (> 1 mas) astrometric microlensing signal is 
detected in OB110462 (§7.1). A filter-dependent astro- 
metric signal is also detected in OB110037 (§7.4); aver- 
aging the astrometry across the F814W and F606W fil- 
ters shows a small but significant 0.4 mas signal. The re- 
maining three targets MB09260 (§7.2), MB10364 (§7.3), 
and OB110310 (87.5) have astrometric signals that are 
below HST’s detection threshold. Table 8 reports the 
lens masses, Einstein radii, maximum astrometric shifts, 
and microlensing parallaxes either as measured values or 
upper limits. 

87.1 - 87.5, the posteriors of the joint photometry and 
astrometry microlensing fits for the individual targets 
are presented. Tables 10 - 13 and Tables 1 - 2 in Lam 
et al. (2022) list the median and lo (68%) credible in- 


tervals for each parameter, in addition to the maximum 
a posteriori (MAP) and maximum likelihood (MLE)'* 
solution. 

Using the methodology described in §6, constraints 
are placed on the lens types of each target. We report 
the relative probabilities of Star:BD:WD:NS:BH in Ta- 
ble 9. Figure 13 shows the distribution of dark vs. lu- 
minous lens probability as a function of mass for each 
target. Table 9 lists the upper limit on the probability 
of a stellar lens and lower limits on the probabilities of 
different dark lenses for each target. In 87.9, we com- 
pare the yield of BHs from our search to that expected 
from simulations assuming there are ~ 10° BHs in the 
Milky Way. 


7.1. OB110462 


The results for OB110462 are presented in the main 
paper, §4.1 of Lam et al. (2022). 


7.2. MB09260 


The data and model for MB09260 are shown in Fig- 
ures 3 (photometry) and 16 (astrometry), and the fit 
posteriors are summarized in Table 10. The inferred 
Einstein crossing time tm is 142.6473? days, the mi- 
crolensing parallax 7 z is Csr, the Einstein radius 
Om is 1.0440-22, and the lens mass is 1.37104 Mo. 

The probability that MB09260 has a dark lens is at 
least 96%, and the probability of a stellar lens is at most 
4%. The relative probabilities of Star:BD:WD:NS:BH 
are 4:0:38:44:14. Stellar lenses are only allowed below 
1Mo. A white dwarf or neutron star is the most proba- 
ble type of lens, with black holes possible but less likely. 


7.3. MB10364 


The data and model for MB10364 are shown in Fig- 
ures 4 (photometry) and 17 (astrometry), and the fit 


posteriors are summarized in Table 11. The inferred 


Einstein crossing time tg is 61.11+};5) days, the mi- 


crolensing parallax 7 is O27 a the Einstein radius 


On is 0.46103}, and the lens mass is 0.217514Mo. 





14 For most of the targets, the MAP solution is equal to the MLE so- 
lution, i.e. values of the parameters for the mode of the posterior 
distribution are the same as those where the likelihood function 
is maximized. 

15 These are the values for the posteriors before they are split into 
modes, which is why the values differ from what is reported in 
Table 10. However, since these parameters are all globally uni- 
modal, their distributions nearly identical across modes within 
the uncertainties, and neither mode is strongly favored, the val- 
ues are not very different and we report the global median and 
uncertainties. 































































































0.35 MB09260 0.35 MB10364 Rae OB110037 
faa «Dark, (41%) Hes «Dark, (9%) Gg: Dark, (0%) 
0.305 SB Darks (55%) 0.305 WBA Darks (55%) 0.30 BBA Darks (26%) 
0.20- 020+ 0-20- 
© 2 © 
a 0.1574 a 0.1574 a 0.1554 
0.10- 0.10- 0.10- 
0.057 0.057 0.05- 
0.00 0.00 0.00 
1 2 3 4 5 0.0 0.2 0.4 0.6 0.8 0.2 0.4 0.6 
Lens mass (Mo) Lens mass (Mo) Lens mass (Mo) 
pa6 OB110310 ge OB110462 DW G6 OB110462 EW 
Gal: «Dark, (18%) Hl «Dark, (91%) G8 Dark, (59%) 
0.305 BR Darks (77%) 0.305 PE Darks (9%) 0.30- BEB Darks (41%) 
0.25- Ha Lum. (5%) 0.25- Ha Lum. (0%) 0.25-4 Ha Lum. (0%) 
0.204 60.204 0.204 
© © © 
a 0.1574 a 0.1574 a 0.1574 
0.107 0.10 0.10- 
0.055 0.057 0.05- 
0.00 - 0.00 0.00 
0 1 2 3 4 3 4 5 1 2 3 4 
Lens mass (Mo) Lens mass (Mo) Lens mass (Mo) 
Figure 13. Probabilities for dark (gray) and luminous lenses (red) , as a function of lens mass for each target. Lenses that 


are dark because all stars of that mass have already evolved and died (implying they are compact objects) are subscripted 
with A (square hatch), while lenses that are dark because a star would be too bright to be hidden in the allowed lens flux are 
subscripted with B (diamond hatch). The probabilities for the luminous lenses are upper limits, while the probabilities for the 
dark lenses are lower limits, since the method described in §6 only places an upper limit on the brightness of the lens. Note 
that there are two fits for OB110462, one with equal weighting to the astrometry and photometry data (OB110462 EW) and 
one with the default weighting of the astrometry and photometry data (OB110462 DW). See §5.3 for details. 


MB10364 is a low mass object, with the possibility 
of a neutron star or BH lens ruled out. The relative 
probabilities of Star:BD:WD are 36:29:35. 


7.4. OB110037 


The data and model for OB110037 are shown in Fig- 
ures 5 (photometry) and 18 (astrometry), the fit poste- 
riors are summarized in Table 12. The inferred Einstein 
crossing time tg is 92.7815°63 days, the microlensing par- 
allax tz is Te is 0.3740°);, the Einstein radius 0g is 
1.241638, and the lens mass is 0.4179:13Mo. 

The probability that OB110037 has a dark lens is at 
least 26%, and the probability of a stellar lens is at most 
74%. The relative probabilities of Star:BD:WD:NS:BH 
are 74:0:26:0:0. Stellar lenses are only allowed below 
0.6Mo, and white dwarfs are the only type of compact 
objects allowed. 

The lightcurve of OB110037 appears to have some 
type of perturbation at MJD ~ 55690. This feature 
is also apparent in the MOA lightcurve, raising our con- 


fidence that the lightcurve feature is real. This pertur- 
bation may be attributed to a binary lens. 

In addition, the astrometry fit, in particular for the 
F606W filter, is quite poor (Figure 18). Although the 
first 5 observations from 2011-2012 seem to agree be- 
tween the two filters, a drastic difference that increases 
as time goes on begins in 2013-2017. This may be at- 
tributed to a binary source. 

It is curious that the photometry may be better ex- 
plained by a binary lens model, while the astrometry is 
likely better explained by a binary source model. Re- 
analysis of this event with both types of binary models 
would be a worthwhile and interesting pursuit, but is 
beyond the current scope of this paper. 


7.5. OB110310 


The data and model for OB110310 are shown in Fig- 
ures 6 (photometry) and 19 (astrometry), and the fit 
posteriors are summarized in Table 13. The inferred 


Einstein crossing time tp is 82.64t7 5" days, the mi- 
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Figure 14. OB110462 astrometry, using the equally weighted likelihood. Left column, top to bottom: RA vs. time with 
maximum likelihood (MLE) unlensed source motion model subtracted; residuals to the MLE model for RA vs. time fit. HST 
F814W astrometry data is shown in red; HST F606W astrometry data is shown in blue. The MLE model is shown in black. 
Fifty random draws from the posterior distribution are shown in light gray. Middle column, top to bottom: Same as left column, 
except Dec instead of RA. Right panel: astrometry as seen on-sky, in the barycentric frame. OB110462 shows a strong > 1 mas 
astrometric microlensing signal. 
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Figure 15. Same as Figure 14, but using the default weighted likelihood for OB110462. 
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Figure 16. MB09260 astrometry. Same as Figure 14, but for MB09260. The astrometric signal is small, and around the limit 
of the precision of the F814W measurements. 
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Figure 17. Same as Figure 14, but for MB10364. The astrometric signal is small, at or below the limit of the precision of the 
F814W measurements. 
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Figure 18. Same as Figure 14, but for 
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OB110037. The photometry and astrometry seem to indicate this object is a binary. 
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Figure 19. Same as Figure 14, but for OB110310. The astrometric signal is small, and around the limit of the precision of the 


F814W measurements. 


Table 10. MB09260 Fit Values, Modes 


















































Parameter Mode 1 Mode 2 
Meat}? MAP MLE Mea}? MAP MLE 
to (MJD) 55099.1971°49 55099.27  55099.27 | 55099.33T] 3S —55099.99 55099.99 
ug —0.0970-09 -0.07 -0.07 0.024803 -0.0002 — -0.0002 
ty (days) 143.1673 33 141.71 141.71 142.3773 33 141.64 141.64 
log19 (8g /mas) 0.037034 0.24 0.24 o.008t eb) 0.06 0.06 
mg (mas) diet 02 0.12 0.12 0.107063 0.09 0.09 
7E,B -0.08t?0n9 -0.08 -0.08 —0.08t0-008 -0.07 -0.07 
7E,N -0.021 0-03 -0.005 -0.005 —0.04 40-08 -0.06 -0.06 
© 50, (mas) 236.2570 13 236.43 236.43 236.2070 13 236.30 236.30 
@§0,N (mas) 692.0775 '49  -692.15 -692.15 | —692.0170-12 — -691.91 -691.91 
Hg, 6 (mas/yr) —5.0070 02 -5.07 -5.07 —4.9970 02 -5.05 -5.05 
Hg, (mas/yr) —3.38 10-08 -3.34 -3.34 —3.3970 02 -3.44 -3.44 
borF.M 0.617803 0.60 0.60 0.607803 0.61 0.61 
™base,M (mag) 17.43 70-003 17.43 17.43 17.43 70-003 17.42 17.42 
OSFF,H8 1,00%:0'02 0.99 0.99 0.997803 1.00 1.00 
™Mbase,Hg (mag) | 17.8440-003 17.83 17.83 17.8470 -003 17.84 17.84 
oS FF,H6 1030-02 1.02 1.02 102.0002 1.02 1.02 
mbase,H6 (mag) | 20.75470-007 20.75 20.75 2007s ey Ope 20.75 20.75 
My (Mg) 1.4470-78 2.70 2.70 T80h0 ies 1.44 1.44 
mp, (mas) oiag.b0-08 0.25 0.25 9:20 9;06 0.20 0.20 
Tel (mas) 0.097 004 0.14 0.14 vig} 05 0.11 0.11 
HL, Rp (mas/yr) —2.6275 38 -0.61 -0.61 —2.887 5°33 -2.82 -2.82 
#1,N (mas/yr) —2.63 70-08 -3.04 -3.04 —2.40+4-09 -1.49 -1.49 
Hret,@ (mas/yr) —2.40+ 0-95 -4.46 -4.46 Sarees -2.23 -2.23 
Hret,N (mas/yr) —0.73 49-89 -0.31 -0.31 —1.0175 38 -1.95 -1.95 
0 (mas) 1.077045 1.73 1.73 ThO2EO Sy 1.15 1.15 
TH 0.09807 0.08 0.08 o.09F 8-04 0.10 0.10 
dc,maa (mas) 0 3a bue 0.61 0.61 0.367014 0.41 0.41 
vw; 0.42 0.58 
logZ 31613.22 31613.55 








NotE—The columns list the median 4 





t1o (68%) credible intervals, maximum a posteriori (MAP) solution, 


and and maximum likelihood estimator (MLE) solution for the microlensing parameters of MB09260. The 
posterior is multimodal (primarily in ug); it has been split and the parameters for each mode reported 
separately. The fraction that each mode contributes to the whole posterior (Sw,) and log evidence 
(log Z) are listed for each mode at the end of the table. 
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Table 11. MB10364 Fit Values 











lo 























Parameter Med? 42 MAP MLE 
to (MJD) 55445.1370°13 55445.06  55445.06 
ug —0.00870-Ot -0.004 -0.004 
tp (days) e111 hoes 61.06 61.06 
log19 (8g /mas) —0.3370-32 -0.40 -0.40 
mg (mas) ai yes 0.11 0.11 
EB —0.24 10-003 -0.24 -0.24 
"BLN 0.12004 0.12 0.12 
© g0,@ (mas) 130.1870 14 130.13 130.13 
«go,N (mas) —78.981 5:14 -79.02 -79.02 
s,m (mas/yr) —7.56 10-06 -7.52 -7.52 
g,N (mas/yr) —6.49t 0-06 -6.47 -6.47 
oSFF,M 0:93 £9°007 0.93 0.93 
™Mbase,M (mag) 15.0219 -00008 15.02 15.02 
oSFF,HS8 0.994803 0.98 0.98 
™base,Hg (mag) 15.327 008 15.32 15.32 
oS FF,H6 1.007803 1.01 1.01 
base, (mag) 16.501 008 16.50 16.50 
My (Mg) gait) 14 0.18 0.18 
mr (mas) daa 0-08 0.22 0.22 
Trey (mas) giz ee 0.11 0.11 
“L,B (mas/yr) 5a4it 8 -5.38 -5.38 
up, (mas/yr) Sarat tee -7.56 -7.56 
rel, (mas/yr) 2467 1:08 -2.13 -2.13 
Mrel,N (mas/yr) ia9t 0-88 1.09 1.09 
@ (mas) 0.464331 0.40 0.40 
TB 0/870 008 0.27 0.27 
5c,maa (mas) 0.167 9-4t 0.14 0.14 








Norr—Same as Table 10, but for MB10364. The solution is uni- 
modal, hence there is only one column and Dw, and log Z are 
not reported. 


crolensing parallax 7 is 0.13193), and Einstein radius 


Om is 0.881095, and the lens mass is 0.781745 Mo.1° 
The probability that OB110310 has a dark lens is at 
least 95%, and the probability of a stellar lens is at most 
5%. The relative probabilities of Star:BD:WD:NS:BH 
are 5:3:65:22:5. Stellar lenses are only allowed below 
1Mo. A white dwarf or neutron star is the most proba- 
ble type of lens, although brown dwarfs and black holes 
are still allowed at the low and high mass ends, respec- 


tively. 


7.6. Source properties inferred from CMDs 


As there is very little blending for all the targets in the 
high resolution F814W and F606W filters (bsrr ~ 1), 
the difference between the target and source on the 
CMD does not change much in color nor magnitude 


16 These are the values for the posteriors before they are split into 
modes, which is why the values differ from what is reported in 
Table 13. However, since these parameters are all globally uni- 
modal, their distributions nearly identical across modes within 
the uncertainties, and neither mode is strongly favored, the val- 
ues are not very different and we report the global median and 
uncertainties. 


Table 12. OB110037 Fit Values 


























Parameter Meat}? MAP MLE 
to (MJD) 55781.5370-25 55781.49  55781.49 
u9 —0.00240-03 -0.008 -0.008 
tp (days) 92.78 13-63 93.31 93.31 
log19 (8 7/mas) 0.0978 th 0.22 0.22 
mg (mas) ovate 2 0.11 0.11 
TELE —0.311 9-008 -0.31 -0.31 
EN agit eer 0.21 0.21 
@g0,@ (mas) 15.211 9-08 15.23 15.23 
@g09,N (mas) —115.5370-02 115.61 -115.61 
wg. (mas/yr) 2.1956 be 2.18 2.18 
g.N (mas/yr) —3.8710-03 -3.86 -3.86 
oSFF.O 0.901008 0.89 0.89 
™base,O (mag) 16.15.+ 0-003 16.15 16.15 
oSFF.H8 0.911008 0.90 0.90 
™Mbase,Hg (mag) 16.33 70-003 16.33 16.33 
oS FF,H6 0.84006 0.82 0.82 
™Mbase,He (Mag) 18.32 70-008 18.31 18.31 
My, (Mo) Gatto 0.55 0.55 
mp (mas) 0.581014 0.74 0.74 
Tet (mas) 0.461014 0.62 0.62 
ep, @ (mas/yr) 6.277} 38 7.59 7.59 
“en, (mas/yr) —6.5670 78 -7.54 -7.54 
rel, B (mas/yr) —4.07+3-37 -5.42 -5.42 
lrel,N (mas/yr) 2.691 0-78 3.68 3.68 
@ 7 (mas) ijaa 0-38 1.67 1.67 
TB a7 308 0.37 0.37 
5c,maa (mas) 0.44 70.13 0.59 0.59 





NotE—Same as Table 11, but for OB110037. 


space (AF814W < 0.1 mag and A(F606W — F814W) 
< 0.1 mag, Figure 2). Based on a CMD analysis, the 
source stars in MB09260, MB10364, OB110037, and 
OB110310 are likely red giant stars in the bulge, as is 
typical for microlensing events in this part of the sky. 

The source of OB110462 in the CMD is around the 
main sequence turnoff (MSTO) on the redder and more 
luminous side of the main sequence, suggesting it is most 
likely a giant or sub-giant star. However, a main se- 
quence source could still be consistent. 

The region of the CMD around the MSTO contains 
both foreground stars from the disk as well as bulge 
stars. We compare the proper motions of OB110462’s 
source to stars in the bright blue foreground as well as 
in the bulge red giant branch to determine which popu- 
lation it most likely belongs to (Figure 20). The source 
is consistent with either population, although it falls 
within the bulk of the bulge population and more on 
the edge of the disk population, hence we consider it is 
most likely a bulge star. This is also consistent with the 
source parallax mg = 0.11 + 0.02 inferred from the fit 
which also indicates the source is likely in the bulge. 

MB09260 and OB110310 are in the the most highly 
reddened field, OB110462 and OB110037 are in inter- 
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Table 13. OB110310 Fit Values, Modes 
































Parameter Mode 1 Mode 2 
Meat}? MAP MLE Meat}? MAP MLE 
to (MJD) 55802.11 74-27 55801.39 55801.39 | 55802.667}°}2 55801.95 55801.95 
uo —0.18 70-08 -0.24 -0.24 —0.00570-08 -0.09 -0.09 
tp (days) 83.401? 39 83.23 83.23 32.20+1-66 81.47 81.47 
log (0 /mas) —0.0570 38 0.04 0.04 —0.06+0-2% 0.04 0.04 
xg (mas) 0.101008 0.10 0.10 6.107 9-02 0.09 0.09 
TELE —0.08t} 01 -0.09 -0.09 [0,097 0s -0.10 -0.10 
EN —0.08t 08 -0.14 -0.14 —o0.1140-h0 -0.21 -0.21 
20, (mas) —104.5670-1% 104.62 -104.62 | —104.6270-09 = 104.58 -104.58 
2g0,N (mas) —183.6170-13 —-183.53 — -183.53 | —183.5770-13 183.49 -183.49 
“g,m (mas/yr) —2.411 9-02 -2.39 -2.39 2.407 0-02 -2.43 -2.43 
g,N (mas/yr) —7.26+ 8-03 -7.26 -7.26 eogt yas -7.28 -7.28 
’SFF,O Ek alias 0.98 0.98 0.961003 0.96 0.96 
™base,O (mag) 18.41 70-008 18.41 18.41 18.4170 -008 18.41 18.41 
OSFF,H8 4.02 3-02 1.04 1.04 1.0270 -02 1.04 1.04 
™Mbase,Hg (mag) | 18.6270-003 18.61 18.61 18.6270-003 18.62 18.62 
oS FF,H6 4027 0-02 1.04 1.04 1.02+0-02 1.05 1.05 
™Mbase, He (mag) | 21.3470-006 21.34 21.34 21.34 40-008 21.34 21.34 
Mz, (Mo) 0.90F 0-77 0.83 0.83 0.71 40-62 0.58 0.58 
mp (mas) Gatto Te 0.28 0.28 0.2370 oF 0.35 0.35 
Tet (mas) Gait eae 0.18 0.18 0.1273 0.25 0.25 
HL, @ (mas/yr) 0.2107 :20 0.22 0.22 —0.16 14-67 -0.30 -0.30 
H1,N (mas/yr) —4.79 42-38 -3.22 -3.22 4672 84 -2.84 -2.84 
Hrel,@ (mas/yr) —2.613-39 -2.60 -2.60 =9.047 7-08 -2.13 -2.13 
lrel,N (mas/yr) —2.4675 98 -4.04 -4.04 ~2.6912-38 -4.44 -4.44 
0 (mas) 0.907848 1.10 1.10 O.a7 7 008 1.10 1.10 
TB oda res 0.16 0.16 Ouiat p08 0.23 0.23 
5c,mazx (mas) 0.321073 0.39 0.39 ogi ae 0.39 0.39 
Dw; 0.43 0.57 
logZ 24631.66 24631.92 
Norr—Same as Table 10, but for OB110310. 
mediately reddened fields, and MB10364 is in the least Osce > Onsr using 
reddened field. This serves as a reminder that within # 
the bulge the amount of extinction is highly variable, bsFF.HST—derived < CUE, (33) 
even over relatively small fields. aki 


7.7. Verifying fit results with bgrr 


The fitting results are validated by comparing the 
best-fit OGLE or MOA I-band source flux fraction 
(bsrr.o, bsrr.m in Tables 1 - 2 in the main paper, and 
Tables 10 - 13) to the high-resolution F814W HST im- 
ages (Figure 12). We only compare the F814W images, 
as it is has a similar effective wavelength to OGLE I- 
band. Assuming a seeing disk radius of ~ 0.65” for 
OGLE and ~ 1.25” for MOA, we add up the flux from 
all stars detected by hst1pass within this radius around 
the target. Next, by assuming all the flux from the 
microlensing target at baseline as seen by HST is due 
to the source (i.e. no flux from the lens or blended 
neighbor stars within HST’s diffraction limited aper- 
ture Ousr ~ 0.09”), we can estimate an upper limit 
on the source flux fraction for a given seeing disk radius 


where F; are the fluxes of the 7 stars within Ose¢ 
and Feyent is the flux of the event at baseline in 
HST. This estimate on the upper limit of the source 
flux fraction derived using HST, bsrro,astT—derived OY 
bsFF.M,HST—derived, is compared to the values inferred 
from the lightcurve fit, bspro or bsrprmu. This ap- 
proach is sufficient to cross check our results, considering 
other differences prevent a perfect comparison (e.g. HST 
F814W is not a perfect match to the OGLE I filter). 

From the photometry fits, MB10364 has bspru ~ 
0.93, OB110037 has bsrro ~ 0.90, OB110310 has 
bsFFO ~ 0.97, and OB110462 has bsFF.O rd 0.05, which 
are all below or no more than ~ 0.05 higher than the up- 
per limits inferred from the HST-derived values shown 
in Figure 12. 

MB09260 has bsprw ~ 0.60, which is higher than 
the upper limit of bs FF,HST—derived ~ 0.4 in F814W 
at 1.25” estimated from the high resolution HST im- 
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Comparison of RGB/Foreground stellar proper motions 
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Figure 20. Top left: CMD with the blue foreground and 
red giant branch (RGB) stars marked. The source color and 
magnitude as inferred from the fit are shown as the yellow 
star. Top right: VPD of the foreground and RGB stars. The 
source proper motion as inferred from the fit is shown as the 
yellow star. Bottom left: Histogram of the foreground and 
RGB RA proper motions. The proper motion of the source is 
shown in the black line. Bottom right: Same as bottom left, 
but for Dec. The source proper motion is consistent with 
either the foreground or RGB bulge population, although 
it falls closer to the bulk of the RGB population, hence we 
consider it is most likely a bulge star. 


ages. However, there is an abrupt transition from 
bsFF,HST—derived ~ 1 to borrusT—derivea ~ 0.4 (in- 
dicating a very bright star) at a radius of 1” from the 
source. The best-fit bsrraz ~ 0.60 falls within these 
two values, indicating that properly treating this situ- 
ation by convolving the HST image with a MOA PSF 
would result in a better match. Hence, MB09260 seems 
consistent with the high resolution image. 

The consistency between the source flux fraction in- 
ferred from the fits and the high resolution imaging pro- 
vides an extra degree of confidence in the inferred fit 
values. 


7.8. Prospects for resolving lens and source 


By definition, microlensing means the source and lens 
are not resolvable. However, after many years, the lens 
and source can separate far enough to be resolved (e.g. 
Batista et al. (2015); Bennett et al. (2015)). In the 
case of dark lenses, such as black holes, non-detections 
of the lens many years after the event can be used to 
place constraints on its properties (e.g. Abdurrahman 


et al. (2021)). For the five targets in this paper, we pro- 
vide estimates to determine if and when taking late-time 
follow-up data would enable such analyses. 

From the results of the fit, we can estimate the time 
necessary to resolve the source and lens tyes via 


tres = Ores / [rel (34) 


where the relative proper motion ju;~1 comes from fitting 
the data, and the minimum angular separation Oe, can 
be estimated using the Rayleigh criterion. For HST with 
a mirror diameter of 2.4m, 0,¢, corresponds to 63.53 mas 
in F606W and 85.35 mas in F814W. 

Assuming there are no contaminating stars, bspr = 
Fs/(Fs + Fy). Then the ratio of lens to source flux, 
or contrast, is F,/Fs = (1 —bsrr)/bsrr. Note that 
bsrr > 1 are allowed by some of the fits which would 
lead to an unphysical negative contrast. bgrr > 1 is of- 
ten referred to as “negative blending” because it means 
there is negative non-source flux (see Equation 25). Neg- 
ative blending can occur in ground-based images if the 
background subtraction is imperfect (Park et al. 2004). 
For HST, where the photometry is not derived from dif- 
ference imaging but rather PSF fitting, this is unlikely 
to be the case. However, when fitting a microlensing 
model to data, negative blending can occur if bs rr © 1, 
simply due to normal photometric uncertainties. This 
is the most likely explanation for HST; all fits with neg- 
ative blending values have posterior probabilities that 
encompass bsrr = 1. Thus, if a contrast value is nega- 
tive, we cap it at 0.'7 In addition, the fit values for pre 
and bgrr are nearly independent, and so these results 
for tres and F',/F's can be thought of as independent. 

The resolving time after photometric peak and max- 
imum contrast in the F814W and F606W HST filters 
for each target are listed in Table 14. Currently, only 
OB110037 could potentially have its source-lens pair re- 
solved in F606W. All the other targets have lens-source 
resolving times at least 15 years post-photometric peak 
or very low contrast. In addition, since they all have 
very high source flux fractions bgrr, their lens-source 
contrasts are very low, which means even after enough 
time has passed for the lens and source to separate, a 
luminous lens would be difficult to detect. Thus, high 
resolution imaging of OB110037 in the near future could 
confirm the results of the fit if the separating source 
and lens could be detected. For the other targets, the 


17 This is technically not correct, analogous to truncating nega- 
tive parallaxes in Gaia; a proper treatment would involve a full 
Bayesian analysis as explained in Luri et al. (2018). Such a treat- 
ment is beyond the scope of this paper. We proceed with capping 
the contrast at 0, as if bs pr ~ 1 then the target is likely very 
dim or dark, and the contrast unlikely to be detectable. 


Table 14. Lens/source resolvability 








Resolving time tres (yr) Contrast (Fz /Fg) 














Parameter F814w F606W F814w F606W 

MB09260 82.08% y py 28.88 pas 0.02903 oo POF 
MB10364 30.721 75°08 oa.87t 38.2? o201 8-02 0.0078 08 
OB110037 17.45 16-82 12.9973-33 0.107807 0.2040 -08 
OB110310 22.037 S555” «16-407 5422) «0009-07. o.004 8-0 
OB110462 Dw i6.91F 7-78 ta.s9t 249 o.o7t8-04 0.034808 
OB110462 EW 15.8168 1i77ts oe 0.0670-0F 9,020.07 











Norre—There are two entries for OB110462, one based on a fit with “default 
weighting” (OB110462 DW) and one based on a fit with “equal weighting” 
(OB110462 EW). Sce §5.3 for details. 


Table 15. Fraction of expected BH de- 
tections vs. tg from PopSyCLE simula- 











tion 
Target tp range % BH 
MB09260 135 < tm < 155 days 50 
MB10364 60 < tp < 62 days 14 
OB110037 87 < tmh< 100 days 12 
OB110310 78 <tpm < 90 days 20 
OB110462 DW 266 < tm < 300 days a7 
OB110462 EW 256 < tm < 325 days 14 





Nots—For each target, the tg range is the median 
+3o0. For OB110462 there are two entries, one 
with equal weighting to the astrometry and pho- 
tometry data (OB110462 EW) and one with the de- 
fault weighting of astrometry and photometry data 
(OB110462 DW). See §5.3 for details. 





absence of any lens detection would imply consistency 
with the results presented here, but could not confirm 
them; however, any detection of a lens would imply the 
fit results here are incorrect. 


7.9. Number of detected BHs 


Next, we compare our observed BH yield to the theo- 
retical expectation calculated using the PopSyCLE simu- 
lations. For a sample of simulated events that would be 
observable by OGLE (see Table 4 of Lam et al. (2020)), 
we calculate the fraction of those events due to BHs as 
a function of the Einstein crossing time tg. Assuming 
the OGLE observability criterion for the MOA sample is 
not strictly correct; however, OGLE and MOA often ob- 
serve an overlapping set of events, so this approximation 
suffices. 

We wish to calculate the probability of detecting k 
BHs in our sample given n events, where the proba- 
bility of detecting a BH in the i-th event is p;. This 
is described by a Poisson binomial distribution, which 
characterizes a “success/no success” experiment with n 
independent trials, where the 7-th trial has probability p; 
of success (Wang 1993). The probability of k successes 
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is given by 


P(K=k)= S° [Iv [[ G-»,) (35) 


AGFRiGA jEAdc 


where Ff, is the set of all subsets of k integers that can 
be selected from the set {1,...,n}, ie. 


Fy, ={A:AC {]1,...,n}, |A] =0,1,..., n}. 


|A| is the number of elements in A, and A® is the com- 
plement of A. In the limit where all p; are equal, the 
Poisson binomial distribution is equivalent to the ordi- 
nary binomial distribution. 

In our case, there are n = 5 independent trials (i.e. 
microlensing targets). We calculate the probability of 
success p; (i.e. BH detection) using PopSyCLE. We define 
the success probability for the 7-th target as the fraction 
of BH lensing events in PopSyCLE over the range of tz 
inferred from the fit 


med(tm) — 30 < tg < med(tg) + 30. 


The probabilities of BH detection for each target are 
listed in Table 15.15 

The results of evaluating Eq. 35 for k = 0,...,5 BH 
detections assuming success probabilities p; is presented 
in Figure 6 in the main paper. PopSyCLE’s Galactic 
model contains 2 x 10° BHs ranging from 5 — 16Mo 
(Lam et al. 2020). There are no 2—5M. mass gap NSs 
or BHs in the simulation, and hence no exact OB110462 
analogue. Thus, we consider a mass gap detection as 
falling between 0 to 1 BH detections. The probability of 
detecting 0 or 1 BHs are ~25% and ~40%, respectively. 
This estimate is consistent with our single detection of 
a NS-BH mass gap object. 


8. DISCUSSION IN DETAIL 


OB110462 is the first definitive detection of a com- 
pact object discovered with astrometric microlensing. 
Depending on the fit, it is either a neutron star (50% 
probability for the EW fit), a black hole (44% probabil- 
ity for the EW fit, 100% probability for the DW fit), ora 
white dwarf (6% probability for the EW fit). MB09260 
and OB110310 are mostly likely white dwarfs or neu- 
tron stars, although a NS-BH mass gap object cannot be 


18 The fact that BH probability does not increase monotonically as 
a function of tg appears to contradict the claims made in the 
main paper that events with longer tg are more likely to have 
a BH lens. The issue is that the uncertainty in BH probability 
rapidly increases at longer tg, as the number of events steeply 
drops off for tg = 50 days (Figure 4 in the main paper). For 
example, for both OB110462 solutions, the estimates are based 
on < 10 events. Because of the small number of events, the 
estimate is highly sensitive to the particular tz range. 
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ruled out. MB10364 and OB110037 are definitively low 
mass objects; OB110037 is most likely a star or white 
dwarf, while MB10364 is either a star, brown dwarf, or 
white dwarf. 

Discussion of the BH yield as compared to theoreti- 
cal expectations is found in §5 of the main paper (Lam 
et al. 2022). Here, we discuss in more detail our sample 
of BH candidates and the questions it raises about the 
Galactic BH population (§8.1), additional observations 
of OB110462 (§8.2), potential sources of systematics in 
OB110462 observations (§8.3), and the future of BH mi- 
crolensing searches (§8.4 and §8.5). 


8.1. OB110462 in comparison to the BH population 


8.1.1. Low-mass X-ray binaries 


Several attempts have been made to determine the 
Milky Way BH mass function using dynamical mass 
measurements of BHs in low-mass X-ray binaries 
(LMXBs). Using a sample of 7 LMXBs, Bailyn et al. 
(1998) found the BH mass function to be tightly cen- 
tered around 7M. with a dearth of systems between 
3 —5Mo; they argued this dearth was not due to ob- 
servational selection effects. Later work by Ozel et al. 
(2010) and Farr et al. (2011) using samples of 15-16 
LMXBs found similar results, providing further support 
to the idea of the NS-BH mass gap. However, Kreidberg 
et al. (2012) cautioned that systematic errors in the anal- 
ysis of LMXB systems could push their inferred masses 
high, artificially creating the mass gap, and Jonker et al. 
(2021) identified potential observational biases that pre- 
vent measurement of high-mass LXMB systems. Addi- 
tionally, LXMBs occupy a very small and specific part of 
BH evolutionary parameter space, and the BHs found in 
those systems are likely not representative of the Galac- 
tic BH population as a whole. 


8.1.2. Filling the NS-BH mass gap 


Measurements of BH masses in non-LMXB systems 
do not show evidence of a mass gap. Gravitational wave 
searches have found mass gap objects both as the merger 
remnant (~ 3M, in GW170917 (Abbott et al. 2017) and 
~ 3.4Me in GW190425 (Abbott et al. 2020a)), and in 
the merger components (~ 2.6Me in GW190814 (Ab- 
bott et al. 2020b)). Non-interacting mass gap BHs of 
~ 3Mo with red giant companions have also been de- 
tected in the Milky Way (Thompson et al. 2019). With 
a mass of 1.6 - 4.4 Mo, the lens of OB110462 is the 
first measured isolated Galactic NS or mass gap BH. 
These detections of mass-gap objects will improve our 
theoretical understanding of compact object formation 
channels. 


8.1.3. A lack of higher-mass BH systems? 


In order to gain a full understanding of the Galac- 
tic BH population, BHs must be uncovered outside of 
closely interacting X-ray binary systems. 

We first consider searches for isolated BHs using mi- 
crolensing, as discussed in this paper. From our sample 
of 5 events, we have a single detection of a mass gap ob- 
ject; all other lenses are lower mass non-BH detections. 
In addition, we have a single detection of a > 1 mas as- 
trometric shift; most of the remaining detections are at 
the ~ 0.5 mas level, near the limit of HST’s precision. 
As discussed in §5.1 of Lam et al. (2022), the low yield 
of BHs in this sample is consistent with predictions by 
Galactic models. However, this presents tentative evi- 
dence that Galactic BHs may be less massive than the 
~ 10Mo expectation. If the BH mass function truly 
peaks at 8Mo, then selecting candidates via long dura- 
tion microlensing events should doubly bias us toward 
finding these high-mass lenses. First, the lensing cross 
section o = 762, is proportional to the mass of the lens 
My, since 0g «x Mz. Thus, more massive objects 
are more likely to be microlenses. Second, the Einstein 
crossing time is proportional to the square root of the 
lens mass, te «x Mz. Thus, long duration events are 
also more likely to be due to massive microlenses. 

Next, we consider searches for BHs in detached/non- 
interacting binary systems. The mass function f of a 
single lined spectroscopic binary is given by 


Popp 3 (1 — €2)3/2 M3 sin? i 
jt Po. G6) 
IG (M, + M2) 





where Pz,» is the orbital period, K is the radial veloc- 
ity (RV) semi-amplitude, e is the orbital eccentricity, i 
is the orbital inclination, MM, is the mass of the visi- 
ble component, and M2 is the mass of the unseen com- 
ponent. If measurements of Pop, K, e, and M, can 
be obtained, then a minimum mass on Mp2 can be de- 
rived. If Mz > 5Me without evidence of luminosity, 
the unseen system is inferred to be a BH. Since sys- 
tems with larger RV semi-amplitudes have larger mass 
functions, they are most likely to host unseen BH com- 
panions. To date, searches for large RV semi-amplitudes 
in spectroscopic catalogs have detected an object which 
fall within the mass gap, but no ~ 10M. BHs, which 
suggest a paucity of higher-mass systems (Thompson 
et al. (2019); El-Badry et al. (2022)). Complementary 
searches using ellipsoidal variables (Rowan et al. 2021) 
also suggest that higher-mass systems are rare. 
Although both microlensing and RV searches should 
be biased toward finding 10M objects more easily than 
3Mo objects, only the latter are being detected. It may 
be that the selection bias for massive objects is cancelled 
out by the fact that the mass function of stars, WDs, 


NS, and BHs sharply decreases from low to high mass. 
Additional work to quantify and compare these two com- 
peting effects, combined with larger sample sizes, will be 
needed to understand the Galactic BH mass function. 

Finally, astrometric searches for detached binaries are 
also eagerly anticipated with Gaia (e.g. Yamaguchi 
et al. (2018); Yalinewich et al. (2018); Wiktorowicz et al. 
(2020)). It will be very fruitful to compare the results of 
those searches to the X-ray transient, microlensing, and 
RV searches. 


8.2. Auxiliary OB110462 observations 


Additional observations of OB110462 may assist in 
ascertaining its NS or BH nature. Here we describe 
planned observations of OB110462, as well as data found 
in searches of archival catalogs. 

The 2021 astrometric observations of OB110462 are 
crucial to the modeling as they extend the temporal 
baseline of the original archival observations by 50%, 
from 6 years to 10 years. The remaining data from HST 
Cycle 29 program GO-16760 to be taken in Fall 2022 
(Lam & Lu 2021) will further extend the baseline by 
another year and improve the characterization of the as- 
trometric signal. The astrometric microlensing shift is a 
deflection with respect to the unlensed position of source 
(Equation 27), and we can only have confidence in our 
measurement of that shift if we also have confidence in 
our measurement of the source’s unlensed proper mo- 
tion. The astrometric shift due to lensing when u > 1 
can be approximated as 6. © Optm/(t — to) (see Ap- 
pendix D). For OB110462 (tg = 280 days and 6g = 4 
mas), the microlensing astrometric shift did not dip be- 
low HST’s astrometric precision of ~ 0.3 mas until 2021, 
10 years after source-lens closest approach. This calcula- 
tion illustrates the importance of having a long temporal 
baseline for OB110462 in order to properly measure the 
source’s unlensed proper motion and characterize the 
astrometric microlensing signal. 

Additional follow-up observations in the X-ray can 
place limits on accretion from the ISM (Agol & 
Kamionkowski 2002). For example, Maeda et al. (2005) 
and Nucita et al. (2006) looked for X-rays at the loca- 
tion of BH microlensing candidate MACHO-96-BLG-5 
reported in Bennett et al. (2002), using ACIS on Chan- 
dra and EPIC on XMM-Newton, respectively. Neither 
detected any X-rays. 

We searched several X-ray catalogs that have observed 
at OB110462’s coordinates to determine whether there 
are any coincident sources. OB110462 was not detected 
as an X-ray source in any of the following catalogs: 
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e Chandra Source Catalog 2.0'9 (Evans et al. 
2019a,b). The limiting sensitivity provides an up- 
per limit of 1.91 x 107 erg/s/cm? at 0.5-7.0 keV. 


e XMM-Newton Science Archive”? (Sarmiento et al. 
2019). This provides an upper limit of 1.52 x 10714 
erg/s/cm? and < 1.38 x 107? counts/s at 0.2 - 12.0 
keV. 


e Swift XRT Point Source Catalogue?! (Evans et al. 
2020). This provides an upper limit of 2.4 x 1073 
counts/s at 0.3 - 10 keV. 


OB110462’s coordinates are not in the eROSITA-DE 
Early Data Release catalog??. 

We also searched the Australia Telescope National Fa- 
cility Pulsar Catalogue”? (Manchester et al. 2005) Ver- 
sion 1.65 for any pulsars coincident with the target. 
There are no coincident pulsars; the nearest pulsar is 
0.55 deg away. 

As mentioned in §3.5, OB110462 is not in the Gaia 
EDR3 catalog. OB114062’s baseline magnitude is 
F606W & 22, while Gaia’s nominal magnitude limit is 
G =& 20.7, and will be even brighter in crowded regions 
(Fabricius et al. 2021). 


8.3. Potential systematics in OB110462 observations 


We briefly discuss whether the tension between the 
photometry and astrometry could be due to systematics 
in the data. 

First we consider the OGLE photometry. Possible 
sources of systematic error include differential refrac- 
tion, proper motion, or low-level stellar variabliliy are 
possibilities for such a long baseline. However, no trends 
due to proper motion or stellar variability are seen in the 
lightcurve. In addition, the F814W — F606W color of 
the bright star within ~ 0.3 arcsec of OB110462 that 
are blended together in the OGLE images are within 
0.1 of each other, and hence effects due to differen- 
tial refraction would be undetectable. We also inspect 
the lightcurves of several stars near OB110462 that 
have similar magnitude to the baseline magnitude of 
OB110462. No trends can be seen that resemble the 
residuals in the EW fit. In addition, we explored re- 
scaling the photometric uncertainties of the OGLE data; 
however, even inflating all uncertainties by a factor of 
~ 2 did not significantly change the structure of the 


19 https: //cxc.harvard.edu/csc/ 

20 http://nxsa.esac.esa.int /nxsa-web/##search 
21 https: //www.swift.ac.uk/2SXPS/ 

22 https: //erosita.mpe.mpg.de/edr/ 


23 https: //www.atnf.csiro.au/research/pulsar/psrcat/ 
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residuals, nor the inferred microlensing parameters, in- 
cluding lens mass. 

Potential sources of systematics in the HST astrome- 
try are discussed in the next section. 


8.4. Improving experimental strategy and design 
8.4.1. Multi-filter astrometry 


In contrast to photometric observations, multi-filter 
astrometry is not routinely obtained, as astrometric ob- 
servations are expensive and facilities with the requisite 
precision are rare. This is one of the first projects to 
explore the impact of different filters on relative astrom- 
etry. The nature of the difference between F606W and 
F814W filter astrometric observations of OB110037 and 
OB110462 is an open question, whether it be astrophys- 
ical (e.g. binaries sources with different color), system- 
atic (e.g. uncorrected CTE), or statistical. However, 
it demonstrates that multi-filter observations are worth 
continuing to pursue in future astrometric microlensing 
studies. For example, it could help break degeneracies 
between certain types of binaries lenses or sources. 


8.4.2. CTE correction 


As mentioned in §4.1, the CTE correction in the _flc 
files is not perfect. Future pursuits will explore other 
methods of correction, such as a re-analysis of OB110462 
that uses the newer and more accurate tabular correc- 
tion for CTE (Anderson 2021). In addition, trying to fix 
CTE via a magnitude-dependent astrometric alignment 
is another avenue that is being explored. 


8.4.3. Observational strategy 


For these precise astrometric measurements, taking 
good observations is critical. The lens mass constraints 
for several of the targets are only upper limits, as the 
astrometric shifts were so small as to be undetectable at 
the precision of the measurements. A dominant source 
of astrometric uncertainty with HST WFC3-UVIS ob- 
servations is the undersampling of the PSF. It has been 
shown that there is a floor in the astrometric precision 
that can be achieved, even at high SNR, when only a few 
exposures are used (Hosek et al. 2015). As the majority 
of the observations in each filter had 4 or less exposures, 
this limited the achievable astrometric precision for sev- 
eral of the targets, in particular OB110462. Increasing 
the number of exposures and implementing small un- 
correlated dithers to sample different pixel phases can 
reduce this floor as VW Ndithers: 

In addition, the effects of CTE are worsening with 
time. Actively mitigating CTE through careful plan- 
ning of observations rather than trying to correct it af- 
terwards is even more important than before (e.g. as 
described in §7 of Anderson et al. (2021)). 


8.4.4. Event selection 


Although all events presented in this work were se- 
lected to have tg > 200 days (Sahu 2009), the inferred 
tz values for 4 of the 5 events did not satisfy this criteria. 
As a result, the true tg range probed extended down to 
tp = 60 days, and did not sample the tg range that max- 
imizes the expected yield of BHs (tg = 100 days). Only 
MB09260 and OB110462 had tg > 100 days, weakening 
the constraints on the BH fraction. We are currently at- 
tempting to determine whether prediction of tg before 
the photometric peak of the event could be improved. 
If possible, this would enable improved selection of BH 
candidates for astrometric follow-up. 

A secondary concern is the target field itself. A suf- 
ficient number of reference stars is needed, hence the 
field in the immediate vicinity of the target must be 
sufficiently crowded. However, the magnitude range of 
those nearby stars must also be similar to the target. Be- 
cause of the steepness of the luminosity function, bright 
targets or targets with high magnification are more dif- 
ficult to analyze as they lack sufficient reference stars to 
perform relative astrometric alignment. This is in ten- 
sion with the need to have high photometric precision 
in order to precisely measure the microlensing parallax. 
Although bright and highly magnified stars should still 
be followed up if they are long duration events, special 
care must be taken when designing observations to en- 
sure good astrometric alignment. 


8.5. BH searches with the Roman Space Telescope 


Although the initial idea and subsequent design re- 
quirements for the Roman Space Telescope (hereafter 
Roman) microlensing survey are driven by exoplanet 
searches (Penny et al. 2019), it also hails the next gen- 
eration of astrometric microlensing campaigns searching 
for BHs. Presently, each event must be followed up indi- 
vidually, with only two facilities (HST and Keck) capa- 
ble of the precision in the near-infrared required to make 
such a measurement. Such measurements are expensive 
(requiring an ~orbit of HST or ~hour of Keck time per 
measurement), prohibiting dense astrometric temporal 
sampling or a large sample of targets. 

Roman will change this with its ability to simultane- 
ously obtain precise photometry and astrometry over a 
wide field of view 100x the area of HST and astromet- 
ric precision almost an order of magnitude better than 
HST (Spergel et al. 2015; WFIRST Astrometry Work- 
ing Group et al. 2019). This will also allow the masses of 
NS-BH mass gap objects to be precisely measured, and 
allow a sample of 100-1000 BH candidates to be built 
up over the duration of the survey (Lam et al. 2020). In 


addition, Roman will probe a large sample of shorter tz 
events which will place constraints on BH kicks. 

To make Roman as effective as possible for finding 
BHs, there are several considerations that need to be 
addressed. Due to the placement of the observatory’s 
solar panels, which dictates the available pointings of 
the telescope, the Galactic Bulge can only be observed 
twice a year during a 72 day window centered on the ver- 
nal and autumnal equinoxes. In addition, for Roman’s 
planned 5 year mission, only 6 of the 10 available Bulge 
seasons are to be dedicated to the microlensing survey. 
These large temporal gaps are generally not a concern 
for exoplanets searches, as the transient portion of the 
lightcurve is nearly covered within the 72 day window. 
However, for long duration events where the transient 
portion of the lightcurve is much longer than the window 
and where a measurement of small microlensing paral- 
lax is crucial, incomplete lightcurve coverage will mean 
the difference between a confirmed BH mass measure- 
ment and only an upper limit. Observations filling in 
these gap will be crucial. Collaboration with a ground- 
based telescope to provide imaging during the gaps, or 
a smaller independent follow-up efforts would be very 
important. 
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APPENDIX 


A. RESCALING OF UNCERTAINTIES 


For each epoch, hstipass returns the RMS error of 
extracted source positions and magnitudes over multiple 


frames oz, Oy and Om, respectively. Following Hosek 
et al. (2015), instead of using the RMS errors for our 
uncertainties, we use the error on the mean a// N where 
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N is the number of frames the star is detected in, inflated 
with a empirical additive error. The empirical constant 
additive error on the positions A,,, and magnitudes A,, 
is calculated for each epoch and added in quadrature to 
the error on the mean. This produces the final rescaled 
uncertainties used in the analysis 


o', = \/o3/N + 3, (A1) 


oy =4/o2/N + 2, (A2) 
Om = VOm/N + Af,. (A3) 


To determine A, and A,,, a sample of bright, unsat- 
urated stars are selected. The exact magnitude range 
constituting “bright” is roughly saturation to 3-5 mag- 
nitudes fainter, with the exact range determined em- 
pirically through the astrometric alignment process. In 
this sample, the additive error is selected such that the 
x? distribution of the reference stars position and mag- 
nitude fits is roughly consistent with the expected dis- 
tribution (Figure 9). The resulting values are listed in 
Table 16. Note that for MB09260 F606W, the same 
additive magnitude error was used across all epochs. 

As mentioned in §4.2.1, for the five microlensing tar- 
gets, the RMS and rescaled astrometric uncertainties are 
generally similar. However, adopting an additive error 
makes the resulting astrometric uncertainties more uni- 
form across the field, particularly for bright stars. This 
is critical as these uncertainties are used as weights in 
the astrometric alignment (§4.2.3 and Figure 9). The 
reference stars’ uncertainties are used to determine how 
good the reference frame transformation is, which is ul- 
timately used to measure the astrometry of the target. 


B. INJECTION AND RECOVERY TESTS 


We use the ks2 software (Anderson et al. (2008); 
Bellini et al. (2018), see also §3 of Sabbi et al. (2016) 
for a detailed description) to inject artificial stars into 
the OB110462 HST images to determine how well we 
can recover the magnitudes and positions of faint sources 
near bright sources. Injection and recovery tests are per- 
formed for all OB110462 epochs in two different man- 
ners: 


1. Sources are injected radially around the bright 
neighbor star at the radius of OB110462 at a 
variety of azimuths excluding the azimuth of 
OB110462 itself, as we cannot recover a source 
planted on top of a real star. 


2. Sources are injected near a star of similar bright- 
ness to the neighbor, in the same radial and 
azimuthal configuration relative to the star as 
OB110462 relative to the bright neighbor. 


Table 16. Additive errors for HST data 








Epoch Pos. error (mas) Mag. error (mmag) 


(yyyy-mm-dd) F814W F606W  F814W F606W 












































MB09260 

2009-10-01 - 0.21 - 4.0 
2009-10-19 0.12 - 3.2 = 
2010-03-22 0.14 - 4.9 — 
2010-06-14 0.14 7 6.1 = 
2010-10-20 0.13 - 9.2 - 
2011-04-19 0.14 0.36 12.8 4.0 
2011-10-24 0.15 0.23 4.0 4.0 
2012-09-25 0.15 0.27 4.3 4.0 
2013-06-17 0.16 0.28 9.8 4.0 
MB10364 

2010-09-13 0.14 1.61 8.6 11.6 
2010-10-26 0.14 0.49 4.5 8.6 
2011-10-31 O17 0.26 5.9 5.0 
2012-09-25 0.18 0.26 8.5 6.1 
2013-10-24 0.17 0.25 6.8 5.8 
OB110037 

2011-08-15 0.11 0.22 16.0 8.5 
2011-09-26 0.11 0.18 4.1 14.0 
2011-11-01 0.10 0.21 4.0 3.8 
2012-05-07 0.10 0.26 12.9 12.2 
2012-09-25 0.10 0.22 6.7 4.3 
2013-10-21 0.10 0.23 3.9 4.9 
2014-10-26 0.10 0.21 4.2 9.7 
2017-03-13 0.10 0.21 4.9 5.1 
2017-09-04 0.11 0.21 8.5 19.6 
OB110310 

2011-09-21 0.10 0.20 9.3 12.9 
2011-10-31 0.12 0.30 4.2 4.9 
2012-04-24 0.13 0.30 10.0 24.3 
2012-09-24 0.14 0.29 10.2 10.9 
2013-10-21 0.13 0.29 4.3 6.3 
2017-03-14 0.13 0.31 4.1 18.8 
2017-09-01 0.14 0.32 10.0 24.9 
OB110462 

2011-08-08 0.13 0.41 18.1 17.0 
2011-10-31 0.19 0.37 3.9 5.4 
2012-09-09 0.19 0.36 4.7 7.2 
2012-09-25 0.19 0.36 6.8 T.2 
2013-10-22 0.19 0.36 3.6 4.8 
2014-10-26 0.19 0.39 3.8 10.7 
2017-08-29 0.20 0.40 5.4 14.7 
2021-10-01 0.17 0.36 7.5 12.9 





The first test directly probes the region around the 
neighbor itself but excludes the actual position of the 
OB110462, while the second test probes a region around 
a star similar to the neighbor, at the same separation 
and angle of OB110462 relative to the neighbor. 

The star we dub “the neighbor” is ~10 pixels (~0.4 
arcsec) west of OB110462. The star we dub “the 
neighbor-like star” in F814W is ~75 pixels (~3 arc- 
sec) northeast of OB110462; for F606W the neighbor- 


like star is ~40 pixels (~1.6 arcsec) north of OB110462 
(Figure 21). 





A6 (arcsec) 














Aa* (arcsec) 


Figure 21. Injection and recovery for OB110462. Black 
circle, center: OB110462. White circle: the bright neigh- 
bor. Blue circle: neighbor-like star for F606W. Cyan circle: 
neighbor-like star for F814W. Positions where stars are in- 
jected in reference to the neighbor-like star are shown as blue 
and cyan stars (F606W and F814W, respectively). 


In summary, both methods produce similar results 
(Figures 22 and 23). We use the single-azimuth results 
from injection around the neighbor-like star in order to 
capture any azimuthal dependence which would be lost 
by averaging over multiple azimuths (Table 17). 

In F814W, the bias is negligible in the first epoch, 
when the target is magnified and is of similar bright- 
ness to the neighbor. However, in later epochs once the 
target is no longer magnified, the measured position of 
the target is biased by ~ 0.4 mas toward the bright 
neighbor star along the target-neighbor separation vec- 
tor. The magnitude and direction of the positional off- 
set are comparable across the two approaches, except in 
2014, where the injection around the neighbor-like star 
leads to a larger difference than when injecting around 
the neighbor. The uncertainties are larger in the ap- 
proach of planting around the neighbor, as it averages 
over more azimuths and results in a larger uncertainty. 

In F606W, the bias is smaller than F814W because 
the shorter wavelength results in higher resolution. Like 
F814W, the bias is also negligible in the first epoch, 
and around ~ 0.25 mas in later epochs. The uncer- 
tainties are also larger when averaging across many az- 
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Table 17. Bias correction derived from injection and re- 
covery 












































Epoch ARA (mas) ADec (mas) AMag (mag) 
F606W 
2011-08-08 -0.025 + 0.001 0.012 + 0.001 -0.001 + 0.000 
2011-10-31 -0.068 + 0.002 0.168 + 0.002 -0.001 + 0.000 
2012-09-09 -0.216 + 0.001 -0.185 + 0.001 -0.005 + 0.000 
2012-09-25 0.108 + 0.002 0.136 + 0.001 -0.003 + 0.000 
2013-10-22 -0.431 + 0.001 0.070 + 0.001 -0.007 + 0.000 
2014-10-26 -0.033 + 0.001 0.298 + 0.001 -0.004 + 0.000 
2017-08-29 -0.144+ 0.001 0.032 + 0.001 -0.011 + 0.000 
2021-10-01 0.124 + 0.002 -0.014 + 0.00 -0.007 + 0.000 
F814W 
2011-08-08 -0.049 + 0.000 -0.014 + 0.000 -0.001 + 0.000 
2011-10-31 -0.203 + 0.001 -0.070 + 0.00 -0.004 + 0.000 
2012-09-09 -0.278 + 0.001 -0.150 + 0.00 -0.016 + 0.000 
2012-09-25 -0.292 + 0.001 0.025 + 0.002 -0.007 + 0.000 
2013-10-22 -0.490 + 0.001 -0.081 + 0.00 -0.011 + 0.000 
2014-10-26 -0.665 + 0.001 -0.138 + 0.00 -0.013 + 0.000 
2017-08-29 -0.258 + 0.001 -0.129 + 0.00 -0.013 + 0.000 
2021-10-01 -0.315 + 0.000  -0.001 + 0.00 -0.007 + 0.000 














NoTEe— Bias correction derived from injection/recovery around a 
star of comparable brightness at the same separation, azimuth, 
and magnitude difference as the target to its bright neighbor. 


imuths. Unlike F814W, the bias is mixed between radial 
and azimuthal components when injecting around the 
neighbor-like star at a single azimuth. 


B.1. Injection around the neighbor 


For each epoch, we determine the separation r of 
OB110462 and the bright neighbor. We inject three 
rings of stars surrounding the bright neighbor, of radii 
r —0.2, r, and r+ 0.2 pixels. Each ring consists of 24 
evenly spaced stars, resulting in one star every 15 de- 
grees. Because we cannot recover stars injected on top 
of the target itself, we do not attempt to recover injected 
stars that fall within 4 pixels of the target itself; this ex- 
cludes three of the positions. We thus inject a total of 
(24 — 3) x 3 = 63 stars per epoch. 


B.2. Injection around the neighbor-like star 


The neighbor-like star we inject around is different 
for F814W and F606W, because the surrounding stars 
do not have the same colors as the neighbor and tar- 
get. The neighbor-like stars were chosen to have sim- 
ilar magnitude and saturation level to the neighbor. 
In F814W, the neighbor tended to be saturated; the 
F814W neighbor-like star is also saturated. On the other 
hand, in F606W the neighbor was not saturated, and 
the F606W neighbor-like star is also not saturated. The 
F814W neighbor-like star is brighter than the neighbor 
in F814W by ~ 0.6 mag, and the F606W neighbor-like 
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star is fainter than the neighbor in F606W by ~ 0.5 
mag. 

In each epoch, we inject three arcs of radii r — 0.2, 
r, and r+ 0.2 pixels centered on the neighbor-like star, 
where r is the target-neighbor separation (~ 10 pixels). 
Each arc consists of 15 stars at the azimuth of the target 
relative to the neighbor + 0.2 pixels/(target - neighbor 
separation in pixels), which corresponds to a subtended 
angle of approximately 2.2 degrees. This corresponds to 
a region of ~0.04 pix? where 3 x 3 x 15 = 135 stars are 
injected. 

At each position we inject stars of magnitude m; — 
0.1, my, and m; + 0.1, where m, is the magnitude that 
results in the same contrast with the neighbor-like star 
as OB110462 to the neighbor. That is, if OB110462 has 
magnitude mr and the neighbor my, and the injected 
star is my; around the neighbor-like star mo, then mr — 





mn = my — MC. 


B.3. Recovery of injected sources 


After planting fake stars into the image, we determine 
how well we can recover the positions and magnitudes. 
To match the properties of our original dataset, we con- 
sider stars to be recovered if they are detected in at 
least N frames, where N is the number of frames that 
were used to calculate the position of the source. Iter- 
ative 3-sigma clipping is performed to exclude outliers 
due to confusion, e.g. from the diffraction spike mask. 
We then use the transformation parameters derived for 
this epoch (as described in §4.2.1) to convert from (x,y) 
pixel positions to (RA, Dec) coordinates. 

We define a polar coordinate system with the origin 
located at the neighbor star for the analysis in §B.1 and 
located at the neighbor-like star for the analysis in §B.2. 
The azimuthal direction is measured counterclockwise 
from the origin-OB110462 separation vector. The aver- 
age offset in the radial and azimuthal directions are r 
and @ (Figures 22 and 23). 

The color of the neighbor and OB110462 are very sim- 
ilar (F606W - F814W = 2.25 and 2.15, respectively). At 
baseline, the neighbor is about 3.1 magnitudes brighter 
than OB110462. During magnification in the first epoch 
(2011-08-08), the neighbor is only about 0.4 magnitudes 
brighter. Since the resolution is higher at shorter wave- 
lengths, it is not unexpected that the positional bias is 
less in F606W than F814W, since the separation is the 
same in both filters. It is also not surprising that the 
bias is smallest in the first epoch when the magnitude 
difference between OB110462 and the neighbor is small, 
and larger in the remaining epochs when the magnitude 
difference is large. 
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Figure 22. Positional bias in F814W (top panel) and 
F606W (bottom panel) as calculated from injection and re- 
covery tests, injecting around a neighbor-like star (left col- 
umn) or around the neighbor (right column). The positional 
bias (position recovered minus the true position input) is 
shown as a function of radial r and azimuthal 0 components 
(top row), as well as total positional bias (bottom row). 


Cc. ABSOLUTE PROPER MOTION REFERENCE 
FRAME 


The HST astrometry in §4.2 is derived in a proce- 
dure where the average motion of the stars is at rest 
with respect to the reference frame. To interpret the 
lens’ proper motion or transverse velocity, we must 
place the relative astrometry into an absolute reference 
frame. To do this, we calculate the proper motion off- 
set between stars in the relative HST frame and the 
absolute Gaia frame. We match all stars in the Gaia 
EDR3 catalog within 1 arcmin of the target that have 
astrometric_excess_noise_sig < 2 and a proper mo- 
tion measurement (i.e. 5-parameter and 6-parameter 
solutions) to the bright stars in our HST proper mo- 
tion catalog (F814W < 22 for MB10364 and OB110037; 
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Figure 23. Magnitude bias in F814W (top panel) and 
F606W (bottom panel) as calculated from injection and re- 
covery tests, injecting around a neighbor-like star (left col- 
umn) or around the neighbor (right column). The magni- 
tude bias is defined as the magnitude recovered minus the 
true magnitude input. 


F814W < 23 for MB09260, OB110310, and OB110462). 
The 1-iteration 30-clipped uncertainty-weighted aver- 
age difference in the proper motion between the cross- 
matched stars is calculated, then applied to the relative 
HST astrometry to place it into the absolute Gaia proper 
motion frame. The values to convert between the HST 
and Gaia frames for each target’s field are listed in Table 
18. 

The vector point diagram of proper motion differences 
between cross-matched sources in the Gaia and HST 
F814W catalogs, after applying the proper motion offset 
to place the HST catalog into the Gaia reference frame, 
is shown in Figure 24. In general, the proper motions of 
bright stars in Gaia are inconsistent with those derived 
using HST. For fainter stars, the uncertainties are much 
larger, so there is more consistency between Gaia and 
HST; however, there is substantial scatter between the 
measurements of the two catalogs. The discrepancies 
between the HST and Gaia proper motions could indi- 
cate that the uncertainties are underestimated in one 
or both catalogs, or that there are higher-order distor- 
tions between the two reference frames that cannot be 
captured by a constant offset. However, the most likely 
explanation is that the Gaia proper motions are not ac- 
curate, as it is clear that the Gaia EDR3 astrometry 


39 


Table 18. HST to Gaia Proper Motion Ref- 
erence Frame Offset 























Target (Ajax) (mas/ys) (Aus) (mas/yr) 
MB09260 -2.56 + 0.13 -4.25 + 0.10 
MB10364 -2.70 + 0.10 -4.56 + 0.10 
OB110037 -2.31 + 0.23 -4.48 + 0.20 
OB110310 -2.20 + 0.12 -4.73 + 0.08 
OB110462 -2.34 + 0.15 -4.72 + 0.11 

NoTe— Aftas = UGaia,ax — HHST F814W,ax and 


Aus = MGaia,s — HST F814W,6- 


is unreliable in crowded regions like the Galactic Bulge 
(see Appendix E and Rybizki et al. (2022)). 


D. ASTROMETRIC LENSING IN GAIA 


As discussed in §4.2.6, there are several differences be- 
tween the proper motion modeling of microlensing tar- 
gets in this paper as compared to Gaia. First, we si- 
multaneously model the proper motions of the lens and 
source with parallax within a microlensing model. In 
contrast, Gaia models the proper motion of a single star 
with parallax. As the source and lens are not resolved by 
Gaia during the observations, if the lens is luminous, the 
proper motion of the source and lens are being treated 
as a single object. Depending on the relative source-lens 
motion and colors, this can produce different inferred 
proper motions. Second, the temporal baseline of the 
data used to calculate the proper motions differ. The 
HST time baseline spans 2010-2013 for MB10364 and 
2011-2017 for OB110037 and OB110310, in contrast to 
2014-2017 in Gaia. As the astrometric lensing signal is 
time dependent, temporal baseline differences can also 
lead to different inferred proper motions. 

With regard to the second point, we estimate the 
amount by which astrometric microlensing could affect 
the proper motion measurement. The lens-source sep- 
aration u(t) (Equation 20) sets the time evolution of 
the astrometric shift 6.(t) (Equation 27). Three terms 
contribute to u(t): the impact parameter uo, the source- 
lens proper motion (t—to)ft,.;/tz, and the microlensing 
parallax te P(t). As the Earth’s orbit around the Sun 
in nearly circular, |P(t)| ~ 1 to within 3%. For typi- 
cal microlensing events like MB10364, OB110037, and 
OB110310, |wo| < 1 and te < 1. Thus, long after 
source-lens closest approach, i.e. (t — to)/tm >> 1, the 
contribution by the impact parameter and microlensing 
parallax terms to u(t) are subdominant to that of the 
source-lens relative proper motion, and the source-lens 
separation can be approximated 


= 


ths 
u(t) & 





> 1. (D4) 
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Figure 24. Vector point diagram of proper motion differences between cross-matched sources in the Gaia and HST F814W 
catalogs within 60” of the targets, after applying the constant proper motion offset derived in Appendix C to place the HST 
catalog into the Gaia reference frame. The color of the point denotes the cross-matched source’s magnitude in Gaia G-band. 
MB10364, OB110037, and OB110310 are sources in Gaia; they are shown as the black-outlined squares and error bars. 


The astrometric shift for large source-lens separations 
u >> 1 can be approximated as 


On 
u(t) 


Putting these two equations together, we can approxi- 
mate how much a proper motion measurement based on 
observations made at times tz and t, will change due to 
astrometric microlensing as 


dc(t2) = dc(t1) 
tg — ty 


0 t t 
~ OF E E . (D7) 
tg—t, \te—to t1 —to 


where (tz — to)/tz, (t1 — to)/tg > 1. For MB10364, 
OB110037, and OB110310, for t; and tg correspond- 
ing to the start and end of Gaia EDR3 observations, 
Apm ®& 0.003, 0.02, and 0.01 mas/yr, respectively. In 
other words, the effect of lensing on proper motions 
should be negligible in Gaia. Thus, if the lens is dark, 
the proper motions of the source as measured by the 
HST data and the proper motion of the target reported 
by Gaia should be the same within the uncertainties. 


de(t) & (D5) 


Apu y (D6) 





E. GAIA DIAGNOSTICS 


We would ideally like to incorporate information from 
Gaia EDR3 into our target analysis (e.g. using the 
reported parallaxes and proper motions to inform the 
prior on the source parallax in the fit, or to compare 
to the posteriors as a cross-check to validate the re- 
sults). However, it is known that there are many as-of- 
yet-unresolved systematics in the Gaia EDR3 astrome- 
try, especially toward the Bulge, and extra care must be 
taken to verify if a proper motion or parallax for a partic- 
ular Gaia source is reliable. Gaia EDR3 is much better 
than Gaia DR2 in terms of photometry, but the astrome- 
try still has issues that need to be worked out (Fabricius 
et al. 2021). We investigate several different metrics for 
the Gaia solutions for our targets to determine the reli- 
ability of the reported parallax and proper motions. A 
brief summary of the meaning of relevant Gaia statistics 
is presented in Table 4; we refer the reader to the Gaia 
Early Data Release 3 documentation (van Leeuwen et al. 
2021) for details. 


E.1. Solution type 


Sources in Gaia have varying amounts 
of information available, described by the 
astrometric_params_solved parameter. A value of 
3 signifies a 2-parameter solution (position), a value of 
31 signifies a 5-parameter solution (position, parallax, 
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Table 19. Distribution of Gaia parameters towards fields 
of interest 











Parameter MB09260 MB10364 OB110037 OB110310 
Nstars 611 1157 1056 884 
fo>2 0.44 0.51 0.59 0.55 
fe=o 0.81 0.77 0.85 0.86 
fmulti=o 0.79 0.41 0.53 0.70 
med(multiso) 4.50 14.00 9.00 5.00 
multi o[95%)) 37.55 54.00 42.30 43.00 
foaa=o 0.76 0.61 0.64 0.64 
med(oddso) 25.50 30.00 26.00 27.50 
odds 995%) ) 83.55 78.20 80.00 78.70 
(log; 9(amp)) -1.24 -0.89 -1.03 -1.13 
Pl0g19(amp) 0.45 0.39 0.58 0.44 
(log, 9 (ruwe)) 0.04 0.12 0.08 0.06 
Plog 19 (ruwe) 0.10 0.16 0.11 0.11 
(logy9(C)) 0.23 0.22 0.25 0.26 
Plogyg(C) 0.12 0.16 0.15 0.15 





NOTE—Nsgtars: Number of Gaia sources within 1 deg? of the target. 
fos2: Fraction of sources with 
astrometric_excess_noise_significance> 2. 
fe=o: Fraction of sources with astrometric_excess_noise> 2. 
fmuitizo: Fraction of stars with ipd_fracmulti_peak = 0. 
med(multis9): Median of nonzero ipd_frac_multi_peak values. 
multisoj95%)): 95th quantile of nonzero ipd_frac_multi_peak values. 
foaa=o: Fraction of stars with ipd_frac_odd_win= 0. 
med(odds0): Median of nonzero ipd_frac_odd_winvalues. 
odds oj95%)): 95th quantile of nonzero ipd_frac_odd_win values. 
fies mie Mean of log, (ipd_gof_harmonic_amplitude). 

Flogy9 (amp)! Standard deviation of 
log;  (ipd_gof _harmonic_amplitude). 

(log, 9(ruwe)): Mean of log; (ruwe). 

Plogyg(ruwe): Standard deviation of log;9(ruwe). 

(log;9(C)) : Mean of log,9(phot_bp_rp_excess_factor). 

Plogy9(C): Standard deviation of log, (phot _bp_rp_excess_factor). 


proper motion), and a value of 95 signifies a 6-parameter 
solution (position, parallax, proper motion, astromet- 
rically estimated effective wavenumber). 5-parameter 
solutions are generally the most accurate, followed by 
6-parameter, then 2-parameter solutions. 6-parameter 
solutions are worse than 5-parameter solutions because 
it means an assumption had to be made about the source 
color. This usually indicates an issue with determining 
the properties of the source, which reduces the accu- 
racy of the solution. This can happen in very crowded 
regions, like the Bulge. MB09260 has a 2-parameter 
solution, MB10364 and OB110310 have 6-parameter 
solutions, and only OB110037 has a 5-parameter solu- 
tion. Based on the high-resolution HST images (Figure 
1) this is not surprising, as OB110037 is the brightest 
object in its vicinity, while the other targets are near 
comparably bright stars. 


E.2. Image parameter determination parameters 


The Gaia EDR3 solution assumes each source 
is a single star. Gaia’s image parameter deter- 
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mination parameters (IPD) can be used to deter- 
mine whether this assumption is valid. Different 
types of binaries can be identified through large 
of IPD values: ipd_gof-harmonic_amplitude (indi- 
cates partially resolved binaries, asymmetric images), 
ipd_frac_multi_peak (indicates resolved, close bi- 
naries), ipd_frac_odd_win (indicates another bright 
source, or observation window conflicts for wide pairs) 
(Lindegren et al. 2021b). Given that the targets all have 
lightcurves well described by point lenses (although see 
§7.4 regarding indications of binarity in OB110037), we 
do not expect any of these parameters to be unusual 
for our targets; unusual values likely indicates problems 
with the astrometric solution itself due to crowding or 
contamination. 

The ipd_gof_harmonic_amplitude values for 
MB09260, MB10364, OB110037, and OB11030 (0.089, 
0.064, 0.036, 0.042, respectively) are somewhat large 
compared to the median values reported in Tables 4-6 
of Lindegren et al. (2021b). However, this is to be ex- 
pected toward the Bulge. Compared to other sources 
within 1 deg? (Table 19) of these targets, the values are 
not unusual.?4 

MB10364 has ipd_frac_multi_peak = 15. The high- 
resolution HST images do not indicate MB10364 is a 
resolved binary. The last HST image is from 2014, and 
Gaia EDR3 observations span 2014-2017. A possibility 
is that the lens and source separation after 2014 was 
large enough to resolve with Gaia EDR3. However, this 
is ruled out by considering the lens-source separation 
criteria (Table 14) and that Gaia’s resolution is worse 
than HST’s. The most likely explanation is confusion 
in Gaia due to stars near MB10364 since the crowding 
within 1 arcsec of MB10364 is high. The other targets 
have ipd_frac_multi_peak = 0. 

MB09260 has ipd_frac_odd_win = 18 and OB110310 
has ipd_frac_odd_win = 55. Similar to the case of 
MB10364 in the previous paragraph, this is likely due 
to bright stars within 1 arcsec of these targets that 
led to confusion in Gaia. The other targets have 
ipd_frac_odd_win = 0. 


E.3. Astrometric goodness-of-fit statistics 


Gaia also presents several relevant goodness- 
of-fit (GOF) statistics from the astrometric fit: 
the renormalized unit weight error (RUWE; ide- 
ally 1), the extra noise required per observa- 
tion to explain the residual in the astrometric 


fit of the source (astrometric_excess_noise; ide- 
ally 0), and the significance of this source noise 
(astrometric_excess noise_sig; insignificant if < 2). 

MB09260 does not have a RUWE because it is not 
calculated for 2-parameter solutions. OB110037 and 
OB110310 both have RUWE ~ 1 (0.971 and 0.981, re- 
spectively). MB10364 has a large RUWE (1.388), al- 
though not unusual compared to the sources nearby (Ta- 
ble 19).?° 

MB09260 and MB10364 both have significant astro- 
metric excess noise (1.241 mas and 0.406 mas, respec- 
tively). OB110037 does not have astrometric excess 
noise and the noise for OB110310 is not significant. 

Table 19 lists the distributions of D = astrometric ex- 
cess noise significance for the stars in Gaia within 1 deg? 
of our targets. Gaia documentation notes that D = 0 
for roughly half the sources and D > 2 for a few per- 
cent of sources with well-behaved astrometric solutions 
(van Leeuwen et al. 2021). Near MB09260, 84% sources 
have D > 0 and 44% have D > 2. Near MB10364, 
79% sources have D > 0 and 51% with have D > 2. 
Near OB110037, 86% sources have D > 0 and 59% have 
D > 2. Near OB110310, 87% sources have D > 0 and 
54% have D > 2. This indicates the astrometric so- 
lution is not well-behaved. This is not surprising, as 
there still exist many systematics in the astrometry, es- 
pecially toward the Bulge. Figure 25 shows astrometric 
excess noise significance plotted against astrometric ex- 
cess noise for systems where astrometric excess noise is 
nonzero. 


E.4. Color excess 


Fabricius et al. (2021) use C = phot_bp_rp_excess_factor 


as a proxy for contamination due to crowding. Stars 
with large excess (C' > 5) tend to have underestimated 
proper motion uncertainties by a factor of ~ 1.7. All of 
our targets have C' < 5; less than 1% of stars within 1 
deg? of the targets have C > 5. This metric does not 
seem to capture crowding toward the Bulge. 


E.5. Summary 


Compared to the high resolution images from HST, 
the various metrics reported by Gaia make sense. 


e OB110037 is the brightest source in its vicinity, 
and hence has a good astrometric solution. Its 
astrometric solution (parallax and proper motion) 
are likely to be reliable. The HST source proper 


24.Note the values reported in Table 19 are the logarithms of 
ipd_gof_harmonic_amplitude. We choose to report the loga- 
rithms as they are approximately normally distributed. 


25 Note the values reported in Table 19 are the logarithms of 
RUWE. We choose to report the logarithms as they are approx- 
imately normally distributed. 


motion agrees with the Gaia proper motion mea- 
surement. 


e MB10364 is bright but is in a crowded field. 
OB110310 is somewhat isolated but is faint. Their 
astrometric solutions are likely to be unreliable. 


e MB09260 does not have enough visibility periods 
to generate an astrometric solution (which requires 
at least 9 visibility periods, while MB09260 only 
has 8). 


In conclusion, although Gaia is a dedicated astromet- 
ric mission, it is not optimized for the crowded and ex- 
tincted Bulge, and the astrometric parameters are likely 
to be untrustworthy there (Fabricius et al. 2021; Rybizki 
et al. 2022). Although Gaia EDR3 is also an improve- 
ment over DR2, those improvements are in the photom- 
etry and not the astrometry. Placing too much weight 
on the Bulge astrometry in analyses (especially the un- 
certainties) should be done with caution. To make use of 
Gaia data for these targets, we will need to wait for fu- 
ture data releases with improved astrometry in crowded 
fields as well as per-epoch astrometry. 


F. GAUSSIAN PROCESS 


The photometric microlensing survey data contains 
uncharacterized systematics in the lightcurves, which 
may be due to unaccounted correlated noise from as- 
trophysical processes or data acquisition and reduction. 
Correlated noise can bias the results of parameter esti- 
mation. Golovich et al. (2022) fit a set of publicly re- 
leased OGLE-III and OGLE-IV microlensing lightcurves 
using gaussian processes (GPs) to account for corre- 
lated noise; we follow that approach here. We use the 
celerite package (Foreman-Mackey et al. 2017) to im- 
plement GPs in our microlensing model. 

For a thorough reference about GPs and their applica- 
tion to machine learning and inference, the reader may 
consult sources such as Rasmussen & Williams (2006). 
In short, a GP is composed of two parts: the noise (the 
stochastic component) and the mean function (the de- 
terministic component). The properties of the stochastic 
component are governed by the covariance matrix, also 
known as the kernel, of the GP. The notation 


y ~ GP(ue(t), Kalt,t)) (F8) 


indicates a GP with mean function 4 with parameters 
6 and kernel K with parameters a. 

When trying to infer some set of parameters 0 from 
data 


y = (y(t) yn (tn))” (F9) 
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the Gaussian log-likelihood is 


1 1 N 
loglL(y|0, a) = —57eKa'Te “5 log |Ka| — 5 log(2z7), 
(F10) 
where N is the number of data points, and rg is the 
residual vector 


yn — pHo(tn))’. (F11) 


In the residual vector, j1g is the mean model, which in 
our case is the microlensing model. The kernel or co- 
variance matrix K, describes the correlated errors. If 
the errors are not correlated, Kq is diagonal and the 
likelihood reduces to the familiar form 


log£(y|@) = -; yy (itor a Gnugto) 


re = (yi — He(t1) 


(F12) 
where o; is the uncertainty on data point y; and y? is 
the “goodness of fit”. 

Following Golovich et al. (2022), a damped driven 
simple harmonic oscillator (DDSHO) kernel kppsHo 
added to a Matérn-3/2 kernel #y73/2 is used to model the 
correlated noise in the photometric microlensing survey 
lightcurves. Both these kernels are stationary, as they 
are a function of the differences of the times only: 


The kernel is given by 
Ka(tij) = KppsHo(Tij) + Km3/2(Tiz) + 5ij0;- (F14) 


Qualitatively, the DDSHO kernel models smooth varia- 
tions, while the Matérn-3/2 captures more irregular vari- 
ations. 

The DDSHO kernel is given by 


—wotis WOTi4 Tv 
KppsHo(Tiz) = Sowoe “° ig /V2 cog (“y = "| 


(F15) 
where So controls the amplitude of the deviation from 
the mean model and wo controls the variation frequency. 
This kernel has been used in asteroseismic modeling (Li 
et al. (2019) and references therein). 

The Matérn-3/2 kernel is given by 


Bray 
KM3/2(Tij) = ge Varii/P (: + Ms) (F16) 


where o determines the amount of deviation allowed 
from the mean model, and p is the characteristic co- 
herence scale. The Matérn-3/2 kernel has been used 
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Figure 25. Astrometric excess noise vs. astrometric excess noise significance as a function of G magnitude, for sources 
with non-zero astrometric excess noise in the 1 deg” fields surrounding the microlensing targets found in Gaia. The tar- 
gets are marked as stars; OB110037 is not marked as it has astrometric excess noise 0. The dotted horizontal line denotes 
astrometric_excess_noise_sig = 2; values < 2 are not significant. 


to model correlated noise in the lightcurves of transit- 
ing and eclipsing exoplanets (Gibson et al. 2013; Evans 
et al. 2015), and in particular is appropriate for modeling 
non-smooth behaviors (Gilbertson et al. 2020). For nu- 
merical reasons (see §4 of Foreman-Mackey et al. 2017), 
the Matérn-3/2 kernel is approximated 


K.M3/2(Tij) =a 





€ 

(F17) 
In the limit € + 0, this is exactly the Matérn-3/2 kernel. 
We implement the approximation with ¢ = 0.01. 

The ground-based OGLE and MOA photometry are 
fit using the Gaussian likelihood with a full covariance 
matrix (Equation F10). The HST photometry and as- 
trometry are fit using the Gaussian likelihood assuming 
a diagonal covariance matrix (Equation F12). 

Note for MB10364, instead of fitting the MOA 
lightcurve using a GP, we instead fit an additive error on 
the ground-based photometry. During nested sampling, 
the GP showed some numerical instability. 


G. PRIORS 


The distributions for the priors 7 are described in this 
sections. N(u,0) denotes a normal distribution with 
mean js and standard deviation 0. Nr(U,0,lo, Uc) de- 
notes a normal distribution with a low end truncation 
at wtol, and a high end truncation at wt+ou,. U(a,b) 
denotes a uniform distribution from a to b. [~+(a, 8) is 
the inverse gamma distribution 


a 


gz *—1 exp|[—8/z]. (G18) 





2 B 
ieee tor = 
(00.8) =F 
The prior distributions for each target are summarized 
in Table 20. 


G.1. Photometry priors 


The five microlensing parameters in a PSPL with par- 
allax fit are to, uo, th, TH n, and Tz,N. 

The prior on to is a normal distribution centered on 
the time of peak magnification in the geocentric frame, 
with a spread of 75 days 


1(to) ~ N (tpeak,@, 75 days). (G19) 


Note that the time at peak magnification in the helio- 
centric frame to is not necessarily the same as in the geo- 
centric frame tpeak,@, hence the large amount of spread 
in the prior. 

The prior on up is a Gaussian with mean 0 and stan- 
dard deviation 0.5 


(uo) ~ N(0, 0.5) (G20) 


(1+*) oe V3801-€) Ti; /p (1-*) evasion . 
€ 
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which takes into account that events with smaller |uo| 
are more likely to be detected, and that events with 
lwo] > 1.5 are not robustly detectable with current 
ground-based surveys.7° 

The prior on tg is a Gaussian centered at 200 days 
with a large spread of o = 100 days. The distribution 
is truncated at —1.80 and 30 (20 and 500 days, respec- 
tively): 


(tz) ~ Np(200, 50, —1.8, 10) days. (G21) 


The priors on the microlensing parallax are estimated 
from the population of bulge microlensing events from 
the PopSyCLE simulation: 


1(tn,B) ~ N(—0.02, 0.12) (G22) 
n(tp,.n) ~ N(—0.03, 0.13). (G23) 


For each dataset filter, bgrr and mpase are also fit. 
For the ground-based photometry, we use a prior 


bsFF ground ~ U(0, 1.1) (G24) 


where the negative blend flux implied by bgrr > 1 al- 
lows for some extra noise such as imperfect background 
subtraction. Similarly for the HST astrometry, we use 
a uniform prior on bs rr 


bsrrusr ~ U(0, 1.05). (G25) 
Mbase 18 a Normal distribution 
Mbase ™ N (Mbase; moses) (G26) 


where Mase is the average magnitude during the un- 
magnified seasons, weighted by the measurement uncer- 
tainties, and om,,,, is 0.1 for OGLE, 0.2 for MOA, and 
0.05 for HST. 


G.2. Gaussian Process hyperparameter priors 


The ground-based photometry includes correlated 
noise we fit. We follow a very similar parametrization 
to Golovich et al. (2022) for the GP priors. The main 
difference is that fit in magnitude space instead of flux 
space, and so our priors are also in magnitudes instead 
of fluxes. 

For oc, we use the prior 


log(a/mag) ~ N(0,5) (G27) 


26 luo| = 1.5 corresponds to a brightening of no more than around 
0.1 mag. When selecting microlensing events, those with a bright- 
ening less than 0.1 mag are generally excluded in survey sam- 
ples to prevent contamination from low-amplitude variables (e.g. 
Mréz et al. (2019)). 
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Table 20. Priors 














Parameter MB09260 MB10364 OB110037 OB110310 OB110462 DW OB110462 EW 
to (MJD) N (55110, 75) N (55460, 75) N (55805, 75) N (55810, 75) N (55770, 75) N (55770, 75) 
uo N (0, 0.5) me z Mm ” 
tr (days) Nr (200, 100, —1.8, 3) ” a » v 
TELE N (—0.02, 0.12) » » » » ” 
TELN N (—0.03, 0.13) » » » » » 
Mbase,O/M (mag) N (17.43, 0.2) N (15.02, 0.2) N (16.15, 0.1) N (18.41, 0.1) N (16.41, 0.1) N (16.41, 0.1) 
bsFF,O/M U(0, 1.1) 2 . x m 2 
Mbase,Hs (mag) N (17.83, 0.05) N (15.32, 0.05) N (16.33, 0.05) N (18.61, 0.05) N (19.85, 0.05) N (19.85, 0.05) 
bsFF,H8 u(0, 1.05) ” " = ” as 
Mbase,H6 (mag) N (20.74, 0.05) N (16.50, 0.05) N (18.33, 0.05) N (21.34, 0.05) N (22.03, 0.05) N (22.03, 0.05) 
bsrr,H6 Uu(0, 1.05) mt ” v ” ” 
log o9,0/m (mag) N (0,5) % 7 " = = 


T~1(0.527, 0.450) - - 
N (1.02e — 03,5) - - 


T~1(0.448, 0.113) 
N(8.41e — 06, 5) 


I~ 1 (0.473, 0.162) 
N (3.60 — 05, 5) 


T~1(0.448, 0.063) 
N(3.53e — 04, 5) 


Poy (days) 

log Wo /mS0,0/M 
(mag? days~?) 

log wo,o/m_(days”*) 
log, 9(@z) (mas) 

mg (mas) 


N (0,5) » » » = 
N (—0.2, 0.3) * » % N (0.5, 0) 
Nr (0.1126, 0.0213, —2.94, 90) * ” » " 
(0.213, 0.250) (0.086, 0.158) U(—0.034,0.091)  U(—0.108,—0.103) — U(0.227, 0.233) U (0.227, 0.233) 
U(—0.697, —0.683) U(—0.096, —0.068) U4(—0.122,—0.104) U(—0.228, —0.154) U(—0.235, —0.183) U(—0.235, —0.183) 
U(—5.96, 1.12) U(—7.93, —1.78) U(0.87, 8.05) U(—1.30, 0.95) U(—4.82, 4.99) U(—4.82, 4.99) 
U(—2.37, 4.17) U(—7.28, 3.41) U(—1.67, 2.94) U(—4.55, —0.49) U(—3.49, 5.91) U(—3.49, 5.91) 





N (0.5, 0) 


@g0,B6 (arcsec) 
&so0,N (arcsec 
Hs,B (mas/yr) 
s,n (mas/yr) 





No te— For definitions of the different variables, see Appendix G. There are two fits for OB110462 depending on the likelihood used, “equal weighting” (OB110462 
EW) or “default weighting” (OB110462 DW). See §5.3 for details. 


which allows a wide range of lightcurve amplitude vari- tm > 50 or te > 120 days as 


ability. 


For p, we use the prior 


7 (logy o(@z/mas)) ~ N(—0.2, 0.3) 
1 (log,9(@2/mas)) ~ N’(0, 0.5), 


(G33) 
(G34) 


p~I7'(a,b) G28) 


where a and b are the constants that satisfy the relation 


med(At) 
0.01 = | [-1(2;a,) dx G29) 
0 





AT 
0.01 = 1 -f [1 (a; a,b) dx G30) 
0 


where med(At) is the median duration between observa- 
tions and AT is the duration of full dataset. This helps 
suppress values at extremely short or long timescales 
that might lead to ill-behaved models.?" 

For So and wo we use the priors 


log Sow) ~ N(med(o?,), 5) (G31) 
log wo ~ N(0, 5). (G32) 


G.3. Astrometry priors 


The prior on the Einstein radius 6g is a lognormal 
distribution estimated from PopSyCLE for events with 


27 See betanalpha.github.io/assets /case_studies/gaussian_processes 


respectively. We use the prior from tg > 120 days for 
OB110462 and tg > 50 days for the other 4 targets. 

The prior on the source parallax 7g is estimated from 
the population of bulge microlensing events from the 
PopSyCLE simulation 


(m5) ~ Np(0.1126, 0.0213, —2.9390, 90.0) mas (G35) 


which corresponds to source distances ranging from 0.5 
to 20 kpc. 

The prior on the source proper motion s,m and LWs,n 
are uniform distributions 


T(uUs,e) ~U(fig p — Sons 5) 4s,p + fozs;,) (G36) 
t(usw) ~Uig n — fons vy Hs,w + fons y) (G87) 


where fig 7, fig,n are the proper motions inferred from 
assuming straight-line motion (no parallax) from the 
F814W data, OTs p> Cs, y are the uncertainties to that 
fit, and f is an inflation factor. To allow a wide range 
of proper motions we use f = 100. 


The prior on the source position at to, Zog,2 and 
ZO0g,N3 is 


T(Xog,z) ~ U(min(rez) — foz,,max(te) + for,) 
(G38) 


1(Xo5,Nn) sod U(min(xy) _ for, ,max(xy) ac FF) 
(G39) 


where rz, Lyn are the positions in the F814W data of 
the target, o2,,, Gx, is the standard deviation, and f is 
an inflation factor. We use f = 5. 


H. ASTROMETRIC COLOR ANALYSIS 


For some stars, the astrometric measurements taken 
in the F814W and F606W filter are discrepant at the 
level of the reported uncertainties. We explore this dis- 
crepancy specifically for OB110462, but this issue is also 
seen in OB110037 and reference stars for all targets. 

First, we must quantify the discrepancy. We consider 
several ways to measure the total offset between the as- 
trometry of the two filters across all epochs. dz is the 
average of the offsets in RA across all epochs for a partic- 
ular star, and can be thought of measuring the amount 
of translation between F814W and F606W. |dz| is the 
average of the magnitude of the offsets in RA across 
all epochs for a particular star, and can be thought of 
measuring the absolute amount of translation between 
F814W and F606W. |dz| is the absolute value of the 
average of the magnitude of the offsets in RA across 
all epochs for a particular star, and can be thought 
of as measuring the total amount of deviation between 
F814W and F606W. Note that |dz| is distinct from |dz]. 
The definitions for dy, |dy], and |dy| are analogous, ex- 
cept they are the offsets in Dec. We also consider these 
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quantities in units of sigma, where the differences in 
each epoch dz; and dy; are normalized by the positional 
uncertainties a, ; and oy ;. 

In Figure 26, we show the distributions of these quan- 
tities as a function of magnitude for stars within 30” of 
OB110462. While not falling in the bulk of the distribu- 
tion, OB110462 is not an extreme outlier. Considering 
how large the variation in positional differences is, es- 
pecially for fainter stars, OB110462 seems well within 
the other positions. For this reason, we assume the 
positional differences are a systematic we can correct 
empirically. We apply a constant positional offset to 
the F606W OB110462 observations to match the posi- 
tions of the F814W observations as described in §4.2.5. 
However, further investigation to determine whether the 
source of the filter dependent astrometry of OB110462 
and other stars may actually be astrophysical is worth 
pursuing, as are more observational programs to study 
filter dependence on astrometry. 

I. DIRECTLY CONFRONTING THE 
PHOTOMETRY AND ASTROMETRY TENSION 


To try and elucidate the tension between the photome- 
try and astrometry, we fit the OGLE photometry alone, 
and separately fit the HST astrometry alone. The re- 
sults are shown in Figure 27. The HST astrometry alone 
does not have very much constraining power-— the uncer- 
tainties on parameters such as to, tg, and uo are so wide 
that the fit itself is not useful. However, the results are 
consistent with those of the EW fit. The OGLE pho- 
tometry has much more constraining power, but alone 
cannot constrain the lens mass. The results are consis- 
tent with those of the DW fit. 
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Figure 26. Difference in astrometry between observations taken in F814W vs. F606W for stars in the OB110462 images, 
aggregated across all epochs. Definitions of dx, |dx|, and |dz| are given in text of Appendix H. The squares mark OB110462. 
The dashed lines indicate the magnitude range of the reference stars in the OB110462 images. 


































0.40 
0.0- == DW Or 
ae 0.30 | (oO een FN 
-02- —- HST ast gf ah ge | 
OGLE phot 7 Oe LON \ 
0.20/ | S . 
\, ; / 
—0.4- 4 C 4 
= 17 =) — 
o E i { / 
< pon shee 
2064 Loe ae ; 
0.1004 oT aee 
0.09 - 
-0.8- 0.08 - a 
\ 7 \ af 
0.07 et. 
0.06 is la 
a4 \ 
T T T T T 0.05 T T T T T T 
0.2 0.0 -0.2 -0.4 -0.6 170 180 190 200 210 220 
rel, E ¢ (deg) 


Figure 27. Left: Comparison of the source-lens relative proper motion inferred from the photometry and astrometry using the 
default weighted likelihood (DW), the photometry and astrometry using the equal weighted likelihood (EW), from the HST 
astrometry alone (HST ast), and from the OGLE photometry alone (OGLE phot). The median is shown as the thick line; the 
shaded region denotes the lo uncertainties. Right: 1 and 20 contours of 7g and direction of the source-lens relative proper 
motion ¢, for the 4 fits presented in the left panel. the EW fit is consistent with the fit to the HST astrometry, although the 
astrometry alone has very little constraining power. The DW and OGLE phot fits are consistent with each other, as expected. 
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